HSUilS 


Mirfloropy  RESOLUTION  TEST  CHART 


I 


I 


AFGUM+09N 


PAPIRA  NO,  MS 


AD-A160  370 

The  AFGL  Global  Spectral  Model:  Expanded 
Resolution  Baseline  Version 


STEPHEN  BftGNNEft 
CHIEN-HSMJNG  YANG 
KENNETH  MITCHEU 


16  November  1984 


Apprpvptf  (m  public 


S 

ft* 


DT1C 

ELECTEj 
OCT  1 7  686 


ATMOSPHERIC  SCIENCES  (XVlSIOff  PROJECT  2310 

AIR  FORCE  GEOPHYSICS  LABORATORY 

HANSCOM  Af^JPA  Pl»| 


85  10 


118 


This  technical  report  has  been  reviewed  and  Is  approved  for  publication' 


FOR  THE  COMMANDER 


DONALD  A.  CHISHOLM,  Chief 
Atmospheric  Prediction  Branch 


SjjHTkt  A.  McCLATCHEY,  Director 
Atmospheric  Sciences  Division 


This  document  has  been  reviewed  by  the  ESD  Public  Affairs  Office  (PA)  and  la 
releasable  to  the  National  Technical  Information  Service  (NTE5). 


Qualified  requestors  nay  obtain  additional  copies  from  the  Defense  Technical 
Information  Center.  All  others  should  apply  to  the  National  Technical 
Information  Service. 


If  your  address  has  changed,  or  if  you  wish  to  be  removed  from  the  mailing 
list,  or  If  the  addressee  Is  no  longer  employed  by  your  organisation,  please 
notify  AF-b/DAA,  Hanscom  AFB,  MA  01731.  This  will  assist  us  In  maintaining 
a  current  mailing  list. 


ia  RtronT  stcufliTv  classification 

Unclassified 


2*  SECURITY  CLASSIFICATION  AUTHORITY 


2b.  DECLASSlFlCATION/DOWNGRAOING  SCHEOULE 


4  PERFORMING  ORGANIZATION  REPORT  NUMBERIS) 

AFGL-TR -84-0308 


6*  NAME  OF  PERFORMING  ORGANIZATION 

Air  Force  Geophysics 


•c.  AOOMS3  iCily.  Sim l*  mnd  ZIP  Codmi 

Hanscom  Air  Force  Base 
Massachusetts  01731 


•p.  NAME  OF  FUNOING/SPONSORING 
ORGANIZATION 


Be  AOORESS  iCily.  Stmt*  and  ZIP  Code) 


REPORT  DOCUMENTATION  PAGE 


lb.  RESTRICTIVE  MARKINGS 


13.  DISTRIBUTIONS  VA I  LABILITY  OF  REPORT 

Approved  for  public  release; 
distribution  unlimited 


6.  MONITORING  ORGANIZATION  REPORT  NUMBERIS) 


7m.  NAME  OF  MONITORING  ORGANIZATION 


7b.  AOORESS  tCity.  Stmt*  mnd  ZIP  Cod*) 


11  TITLE  i Include  Security  Clateificatiom 

The  AFGL  Global  Spectral  Model:  (Contd) 


12.  PERSONAL  AUTMORIS) 

Brenner,  Stephen;  Yang,  Chien-Hsiung;  Mitchell*  Kenneth 


9.  PROCUREMENT  INSTRUMENT  IDENTIFICATION  NUMBER 

1 10  SOURCE  OF  FUNDING  NOS 

PROGRAM 

PROJECT 

TASK 

WORK  UNIT 

ELEMENT  NO 

NO 

NO 

NO 

62101F 

2310 

G7 

03 

13b  TYPE  OF  REPORT 

Final 


IB.  SUFPLEMENTARV  NOTATION 


13b  TIME  COVERED  I  14  DATE  OF  REPORT  (Yr  .  Mo..  Doyl 

from  1/10/82  to  30/9/84  1984  November  16 


COSAT l COOES 


18  SUBJECT  TERMS  fConfinu#  on  rwi*r«#  ,f  nectatary  and  identify  by  block  number) 

Atmospheric  models  Numerical  weather  prediction  , 

Spectral  models  Global  atmospheric  circulation 


19.  ABSTRACT  rConfmu*  on  revert*  if  neettsmry  and  identify  6y  block  number i 

-  N  The  flexible  resolution/truncation  baseline  version  of  the  AFGL  global  spectral 
model  as  adapted  to  the  CRAY-1  is  described.  A  series  of  low -re solution  (6  layer, 
rhomboidal  15)  and  high -resolution  (12  layer,  rhomboidal  30)  forecasts  were  run  and 
compared  to  test  the  performance  of  the  model.  In  general,  higher  resolution  resulted 
in  improved  forecast  skill  in  the  24-to-96-hour  range.  The  only  exception  to  this  is  the 
humidity  forecast,  which  shows  minimal  skill.  This  characteristic  is  rather  insensitiv* 
to  the  resolution  partly  because  of  the  poor  quality  of  analyzed  humidity  fields  used  for 
initial  data  and  verification.  .  ■i-') 

The  original  gridded  (2.  5*  x  2.  5”)  topography  has  been  replaced  by  a  smoothed  ter¬ 
rain  field  that  has  been  passed  through  a  nine -point  smoother,  interpolated  to  the 
model's  Gaussian  grid,  and  then  spectrally  truncated. 

Finally,  the  effects  of  initialization  have  been  studied  by  comparing  a  series  of 
forecasts  subjected  to  several  initialization  methods.  For  forecasts  beyond  (Contd) 


20  DISTRIBUTION/AVAILABILITY  OF  ABSTRACT 
UNCLASSIFIED/UNLIMITED  S  same  AS  RPT  c  DT1C  USERS  □ 


21  ABSTRACT  SECURITY  CLASSIFICATION 

Unclassified 


22*.  NAME  OF  RESPONSIBLE  INDIVIDUAL 

22b  TELEPHONE  NUMBER 
i Include  Area  Cod*i 

22c  OFFICE  SYMBOL 

Stephen  Brenner 

(617)  861-2962 

LYP 

DO  FORM  1473,  83  APR 


EDITION  OF  1  JAN  73  iS  OBSOLETE 


Unclassified _ _ 

SECURITY  CLASSIFICATION  OF  THIS  PAGE 


%  V 


Block  11  (Contd): 

Expanded  Resolution  Baseline  Version 
Block  19  (Contd): 

24  hours,  the  model  is  able  to  supress  spurious  gravity  waves  through  the  combined 
effects  of  the  semi-implicit  time  scheme  and  the  subgrid  scale  diffusion.  The  impact  of 
normal  mode  initialization  is  seen  mainly  in  the  very  short-range  forecasts  (less  than 
24  hours)  and  is  thus  important  for  providing  smooth  first-guess  fields  for  the  analysis/ * 
data  assimilation  cycle. 


msmnisj 


Contents 


1.  INTRODUCTION  1 

2.  EXPANDED  RESOLUTION  BASELINE  MODEL  1 

3.  FLEXIBLE  RESOLUTION/TRUNCATION  4 

4.  VECTORIZATION  8 

5.  REVISED  PREPROCESSING  AND  POSTPROCESSING  10 

6.  MODEL  PERFORMANCE  -  COMPARISON  OF  HIGH-  AND 

LOW -RESOLUTION  FORECASTS  12 

7.  REPRESENTATION  OF  TOPOGRAPHY  38 

8.  EFFECTS  OF  INITIALIZATION  55 

9.  CONCLUSION  64 

REFERENCES  69 

LIST  OF  SYMBOLS  71 


Illustrations 


1.  Spectral  Domain  Showing  Low  and  High  Resolutions,  Rhomboidal 

and  T  riangular  T  runcations  6 


'*.*  "Ti*  rrf*_J 


^T.  1 TJirjirvrj! 


ell 


> 


Illustrations 


2.  Arrays  for  Storing  Spectral  Coefficients:  (a)  Rhomboidal 

Truncation*  and  (b)  Triangular  Truncation 

3.  Schematic  Diagram  of  Model  Evaluation 

4.  Meridional  Profile  of  Processing  Error  -  RMS  Error  of  500  mb 

Heights  for  Various  Resolutions  (All  Rhomboidal  Truncation) 
for  FGGE  Data  of  00Z  18  January  1978 

5.  Meridional  Profile  of  Processing  Error  -  RMS  Error  of  250 mb 

Temperatures  for  Various  Resolutions  (All  Rhomboidal 
Truncation)  for  FGGE  Data  of  00Z  18  January  1978 

6.  Growth  Rate  of  the  Average  1000-200  mb  Height  Errors  for 

Six  Forecasts  (Three  From  January  1978  and  Three  From 
July  1978) 

7.  Growth  Rate  of  the  Average  850  mb  Relative  Humidity  Errors 

for  Six  Forecasts  (Three  From  January  1978  and  Three  From 
July  1978) 

8a.  48-Hour  Forecast  of  1000  mb  Heights  Beginning  From  00Z 

17  January  1978:  Verifying  Analysis.  Contour  interval  is  50  m 

8b.  48-Hour  Low -Resolution  (15R6)  Forecast  of  1000  mb  Heights 

Beginning  From  00Z  17  January  1978.  Contour  interval  is  50  m 

8c.  48 -Hour  High-Resolution  (30R12)  Forecast  of  1000  mb  Heights 

Beginning  From  00Z  17  January  1978.  Contour  interval  is  50  m 

8d.  48-Hour  Forecast  of  1000  mb  Heights  Beginning  From  00Z 

17  January  1978:  Difference  Between  High  and  Low  Resolutions. 
Contour  interval  is  50  m 

9a.  48  -Hour  Forecast  of  500  mb  Heights  Beginning  From  00Z 

17  January  1978:  Verifying  Analysis.  Contour  interval  is  50  m 

9b.  48-Hour  Low-Resolution  (15R6)  Forecast  of  500  mb  Heights 

Beginning  From  00Z  17  January  1978.  Contour  interval  is  50  m 

9c.  48 -Hour  High -Resolution  (30R12)  Forecast  of  500  mb  Heights 

Beginning  From  00Z  17  January  1978.  Contour  interval  is  50  m 

9d.  48 -Hour  Forecast  of  500  mb  Heights  Beginning  From  00Z 

17  January  1978:  Difference  Between  High  and  Low  Resolutions. 
Contour  interval  is  50  m 

10a.  96-Hour  Forecast  of  1000  mb  Heights  Beginning  From  00Z 

17  January  1978:  Verifying  Analysis.  Contour  interval  is  50  m 

10b.  96-Hour  Low-Resolution  (15R6)  Forecast  of  1000  mb  Heights 

Beginning  From  00Z  17  January  1978.  Contour  interval  is  50  m 

10c.  96-Hour  High -Resolution  (30R12)  Forecast  of  1000  mb  Heights 

Beginning  From  00Z  17  January  1978.  Contour  interval  is  50  m 

lOd.  96-Hour  Forecast  of  1000  mb  Heights  Beginning  From  00Z 

17  January  1978:  Difference  Between  High  and  Low  Resolutions. 
Contour  interval  is  50  m 

11a.  96-Hour  Forecast  of  500  mb  Heights  Beginning  From  00Z 

17  January  1978:  Verifying  Arwdysis.  Contour  interval  is  50  m 

lib.  96-Hour  Low-Resolution  (15R6)  Forecast  of  500  mb  Heights 

Beginning  From  00Z  17  January  1978.  Contour  interval  is  50  m 


7 

17 

18 


19 

20 


21 

22 

23 

24 

25 

26 

27 

28 


29 

30 

31 

32 


33 

34 

35 


. 


iv 


Illustrations 


11c.  96-Hour  High-Resolution  (30R12)  Forecast  of  500  mb  Heights 

Beginning  From  00Z  17  January  1978.  Contour  interval  is  50  m  36 

lid.  96-Hour  Forecast  of  500 mb  Heights  Beginning  From  00Z 

17  January  1978:  Difference  Between  High  and  Low  Resolutions. 

Contour  interval  is  50  m  37 

12.  Growth  Rate  of  Height  Errors  for  Forecasts  From  00Z 

17  January  1978:  (a)  1000  mb,  (b)  500  mb,  (c)  250  mb,  and 

(d)  Average  1000  -  200  mb  39 

13a.  Definitions  of  Quantities,  Processes,  and  Errors  in  Representation 

of  Topography  Using  Trapezoidal  Quadrature  42 

13b.  Definitions  of  Quantities,  Processes,  and  Errors  in  Representation 
of  Topography  Using  Gaussian  Quadrature,  Data  on  Equally 
Spaced  Latitudes  43 

13c.  Definitions  of  Quantities,  Processes,  and  Errors  in  Representation 
of  Topography  Using  Gaussian  Quadrature,  Data  on  Gaussian 
Latitudes  44 

14a.  Meridional  Profiles  of  Terrain  Parameters:  DM-Zonal  Mean, 

DE-Zonal  Mean  Plus  Standard  Deviation,  DMAX-Maximum 

Value,  for  Unsmoothed  Terrain  46 

14b.  Meridional  Profiles  of  Terrain  Parameters:  DM-Zonal  Mean, 

DE-Zonal  Mean  Plus  Standard  Deviation,  DMAX-Maximum 

Value,  for  Smoothed  Terrain  46 

15a.  Vertical  Profiles  of  the  Group  Means  of  the  Global  RMS  of  the 
Processing  Error  at  Day  0  and  the  Forecast  Errors  at 

Days  1,  2,  and  3,  With  Unsmoothed  Topography  53 

15b.  Vertical  Profiles  of  the  Group  Means  of  the  Global  RMS  of  the 
Processing  Error  at  Day  0  and  the  Forecast  Errors  at 
Days  1,  2,  and  3,  With  Smoothed  Topography  53 

16.  Vertical  Profiles  of  the  Group  Means  and  Standard  Deviations 

of  the  Differences  Between  the  Global  RMS  Errors  Using 
Different  Terrain  Fields:  (a)  Processing  Error  at  Day  0, 
and  Forecast  Errors  at  (b)  Day  1,  (c)  Day  2,  and  (d)  Day  3  54 

17.  Growth  Rate  of  Vertically  Integrated  Global  RMS  Divergence 

for  High-Resolution  Forecasts  Using  Different  Methods  of 
Initialization,  Beginning  From  00Z  1 7  January  1978 

(FGGE  III -A,  Hough  Analysis)  58 

18.  Vertical  Profiles  (Model  o -Layers)  of  RMS  Divergence  for 

High-Resolution  Forecasts  Using  Different  Methods  of 
Initialization,  Beginning  From  00Z  17  January  1978: 

(a)  Initial  Values,  (b)  12  Hours,  (c)  24  Hours,  and 

(d)  48  Hours  59 

19.  Growth  Rate  of  Vertically  Integrated  Global  RMS  Divergence 

for  High-Resolution  Forecasts  Using  Different  Methods  of 
Initialization,  Beginning  From  12Z  17  February  1979 

(FGGE  III -A,  Optimal  Interpolation)  60 


v 


Illustrations 


20.  Vertical  Profiles  (Model  <7-Layers)  of  RMS  Divergence  for 

High-Resolution  Forecasts  Using  Different  Methods  of 
Initialization:  (a)  Initial  Values,  and  Forecast  Errors 
Beginning  From  12Z  17  February  1979,  (b)  2  Hours, 

(c)  6  Hours,  (d)  12  Hours  62-63 

21.  Vertical  Profiles  of  RMS  Divergence  for  High -Resolution  Data 

on  12Z  18  February  1979:  (a)  Mandatory  Pressure  Levels, 

and  (b)  Model  a  Layers  64 

22a.  Difference  Between  Initialized  and  Uninitialized  Velocity  Potential 
at  <r=  .5  for  Forecasts  Beginning  From  12Z  17  February  1979: 

6-Hour  Forecast.  Contour  Interval  is  5  x  105  m1 2 3 4 5 6 7 * 9s'l  66 

22b.  Difference  Between  Initialized  and  Uninitialized  Velocity  Potential 
at  a=  .5  for  Forecasts  Beginning  From  12Z  17  February  1979: 

48 -Hour  Forecast.  Contour  Interval  is  5  x  10^  m2s"l  67 


Tables 

1.  Computational  Requirements  on  the  CRAY-1  for  Various 

Resolutions  of  the  Moist  Model  9 

2.  Effects  of  Restructuring  Code  (Vectorization)  on  Execution 

Time  for  a  24-Hour  Moist  Forecast  Using  Six  Layers  and 

Rhomboidal  M  =  15  9 

3.  Pressure  Values  (mb)  of  Levels  to  Which  Layer  Temperatures 

Are  Assigned.  The  NMC  12-layer  a  structure  is  used  13 

4.  The  a -Structures  of  the  Low  -  and  High-Resolution  Models  14 

5.  Performance  (Global  RMS  Differences  and  Standard  Deviations), 

January  1978  15 

6.  Performance  (Global  RMS  Differences  and  Standard  Deviations), 

July  1978  16 

7.  Amounts  of  Variance  in  Various  Spectral  Ranges  in  the  Unsmoothed 

and  Smoothed  Terrains  (Units  of  m2)  45 

8a.  RMS  of  the  Error  in  Synthesis  (E^)  and  Error  in  Reproduction  (E2) 

With  Different  Rhomboidal  Truncations  of  Unsmoothed  Terrain 

(Units  of  m)  48 

8b.  RMS  of  the  Error  in  Synthesis  (Ej)  and  Error  in  Reproduction  (E2) 

With  Different  Rhomboidal  Truncations  of  Smoothed  Terrain 

(Units  of  m)  48 

9.  RMS  of  the  Error  in  Synthesis  (E'i)  and  Error  in  Reproduction  (E'g) 

With  Different  Rhomboidal  Truncations  in  Reference  to  the  Given 
Values  on  the  Gaussian  Latitudes  of  the  Unsmoothed  and  Smoothed 
Terrain  (Units  of  m) 


49 


10a.  The  Square  Root  of  the  Power  of  the  Error  in  Reproduction  in 
the  Spectral  Domain  (E3  or  E’3)  for  the  Unsmoothed  Terrain 
(Units  of  m) 

10b.  The  Square  Root  of  the  Power  of  the  Error  in  Reproduction  in 
the  Spectral  Domain  (E3  or  E'3>  for  the  Smoothed  Terrain 
(Units  of  m) 

11.  Definition  of  Spectral  Ranges  Used  in  Tables  12a  and  12b. 
Given  are  the  lower  and  upper  bounds  for  each  index,  K 

12a.  Total  Power  (in  m^)  and  Spectral  Composition  (in  Percent)  of 
E3  for  Trapezoidal  Quadrature 

12b.  Total  Power  (in  m^)  and  Spectral  Composition  (in  Percent)  of 
Eg  for  Gauss-Legendre  Quadrature 

13.  Summary  and  Acronyms  of  Different  Combinations  of  Analysis 
and  Initialization  That  Were  Tested 


The  AFGL  Global  Spectral  Model:  Expanded 
Resolution  Baseline  Version 


1.  INTRODUCTION 

This  report  is  the  second  in  a  series  that  will  serve  as  the  documentation  of 
the  AFGL  global  spectral  model*  The  progress  reported  here  deals  primarily  with 
(1)  the  changes  necessary  to  expand  the  resolution  of  the  baseline  model1  and  adapt 
it  to  the  CRAY-1  computer,  and  (2)  evaluation  of  the  high -resolution  (30R12)  model. 
The  term  "baseline  version"  version  refers  to  the  mathematical /numerical  formu¬ 
lation  and  model  configuration  as  described  by  Brenner  et  al 1  (hereafter  referred 
to  as  TR-1).  The  reader  is  referred  to  that  ’■eport  for  the  appropriate  details. 

2.  EXPANDED  RESOLUTION  BASELINE  MODEL 

As  mentioned  in  section  1,  the  expanded  resolution  version  of  the  AFGL  global 
model  is  an  extension  of  the  baseline  version  described  in  TR-1.  Briefly,  the 
model  is  a  global  spectral  model  using  spherical  harmonics  to  represent  horizontal 
variations  and  discrete  layers  with  the  Arakawa  differencing  scheme  to  represent 
vertical  variations.  Nonlinear  terms  in  the  tendencies  of  vorticity,  divergence, 

(Received  for  publication  15  November  1984) 

1.  Brenner,  S.,  Yang,  C.,  and  Yee,  S.  (1982)  The  AFGL  Spectral  Model  of  the 
Moist  Global  Atmosphere:  Documentation  of  the  Baseline  Version, 

AFGL -TR -82 -03 93,  Ad  A129283. 


<■.  v.  ■*-. 1 


■*,  •*.  ■*.  .  •*.  • 


temperature,  moisture,  and  surface  pressure  are  computed  using  the  transform 

2  3 

method.  '  Time  differencing  is  performed  with  the  semi-implicit  scheme.  The 
independent  variables  of  the  model  are  time,  t;  longitude  and  latitude.  A  and  *;  and 
the  terrain  following  vertical  coordinate,  a=  p/p*.  Dependent  variables  of  the 
model  include: 

(1)  prognostic  variables 
T)  absolute  vorticity 
D  divergence 

T  temperature 

q  natural  log  of  surface  pressure 
Q  specific  humidity 

(2)  diagnostic  variables 
*  geopotential 

a  vertical  velocity 

U,  V  horizontal  pseudo -velocity  components 
For  completeness,  the  continuous  model  equations  are  given  here: 
the  vorticity  and  divergence  equations, 

dV  1  fc>A  ,  „ 

—  - 5—  — —  +  cos *  — - -  +  k'V  x  F,  (1) 

a  J  ^  ^ 

—  =  - - - [  22.  -  cos*  -  V2  Ie  +  *  +  RT0ql  +  V  •  Fh;  (2) 

dt  a  cos2*  IdA  d*J 


the  thermodynamic  equation. 


___L_  auri+  c0s*^TL  +  T,D 

a  cos^*  dk  o<p 


o'i  |I2-  t«Ifc-C-5l  +  JL+hFT: 
da  cp  l  1  cp  h  T 


2.  Eliasen,  E.,  Machenauer,  B.,  and  Rasmussen,  E.  (1970) 

Method  for  Integration  of  the  Hydrodynamical  Equations  With  a  Spectral  Re 


_ onzontai  rieias,  insi.  01  1  neoret. 

Copenhagen.  Report  No.  2. 

Orszag,  S.A.  (1970)  Transform  method  for  calculation  of  vector  coupled  sums: 
Application  to  the  spectral  form  of  the  vorticity  equation,  J.  Atmos.  Sci. 
27:890-895. 


the  moisture  equation. 


IQ.  =  . _ *_[£^Q  +  QD  -<r|2.+  S  +  hFQ; 

dt  a  cos 2</>  dk  d<t>  da  n  w 


the  continuity  equation, 


.2fL=  -  C  -  D; 

at 


the  hydrostatic  equation. 


PQ  =  -  c  T<7  *' 
da«  P  * 


and  diagnostic  equations  for  the  velocity  components. 


-a  -a 

a  =  a  (C  +  D)  -  C  -  D 


n  -  .  cos#  a# a.  j_  ax 

a  e<A  a  c>X 


i  a^+  cos^>  ax 

a  dX  a 


f  =  V2* 

D  =  V2y 


where 


A,TOt«,|v,arco^ffltipL 

da  a  d<p  p*  oa 


aa  a  ax  p*  a<7 


U2  +  V2 
2  cos 


.N 


q  *  lnp* 

(15) 

C  «  v'Fq. 

(16) 

-a  . 

(  )  ■  JQ(  )  doand  (~)  - 

1 

0(  )  do. 

(17) 

Fh  is  the  parameterized  subgrid-scale  horizontal  diffusion*  and  T'  is  the  deviation 
of  temperature  from  the  reference  value  Tq(o). 

Parameterized  physical  processes  in  the  model  include  the  following:  a  bulk 
aerodynamic  boundary  layer  formulation,  large-scale  condensation.  Kuo-type 
moist  convection,  dry  convective  adjustment,  and  subgrid  scale  diffusion.  Once 
again,  the  reader  is  referred  to  TR-1  for  the  details  of  these  aspects  of  the  model. 


3.  FLEXIBLE  RESOLUTION/TRUNCATION 

The  baseline  version  of  the  AFGL  global  spectral  model  was  originally  coded 
for  the  CDC-6600  with  a  user  available  central  memory  of  approximately  98  kilo- 
words  (K-words).  To  keep  the  model  core  contained,  it  was  configured  with  six 
layers  and  rhomboidal  truncation  at  wavenumber  15,  thus  requiring  85  K-words  for 
execution. 

At  the  time,  we  realized  that  the  baseline  version  was  only  a  preliminary  step 
in  the  construction  of  the  model.  For  future  research,  especially  in  the  area  of 
cloud  forecasting,  a  higher  resolution  model  would  be  needed,  one  that  was  at 
least  comparable  to  the  current  operational  spectral  model  at  the  National  Meteor- 

4 

ological  Center  (NMC)  (that  is,  12  layers  and  rhomboidal  truncation  at  wavenum¬ 
ber  30). 

To  construct  the  high -resolution  model,  the  obvious  procedure  would  have 
been  to  replace  all  of  the  resolution -dependent  parameters  (for  example,  array 
dimensions  and  loop  limits)  with  appropriate  high-resolution  values.  However, 
this  again  would  limit  the  model  to  a  particular  resolution.  Thus,  to  make  the 

4 

model  as  flexible  as  possible,  following  the  approach  of  Sela,  the  resolution-de¬ 
pendent  parameters  were  replaced  with  dummy  non-Fortran  symbols  that  could  be 
easily  changed  to  any  values  using  a  text  editor  or  some  type  of  Fortran  prepro¬ 
cessor.  The  same  technique  was  also  applied  to  the  data  preprocessing,  forecast 
postprocessing,  and  normal  mode  initialization  codes. 

In  coding  the  model,  two  additional  aspects  of  flexibility  were  also  considered, 

4.  Sela,  J.  (1980)  Spectral  modeling  at  the  National  Meteorological  Center,  Mon. 
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modular  structure  and  portability  of  the  code*  In  the  area  of  modularity,  as  much 
as  possible,  the  code  was  structured  into  subroutines  corresponding  to  various 
mathematical  and  physical  aspects  of  the  model.  This  allows  for  easy  replacement 
of  various  components  of  the  model  with  new  algorithms.  This  flexibility  is  cru¬ 
cial.  especially  in  the  parts  of  the  model  that  handle  the  physical  parameterizatlons 
that  will  be  replaced  with  new  formulations.  To  keep  the  model  as  portable  as 
possible,  all  machine -dependent  structures  and  subroutine  calls  were  replaced  by 
standard  Fortran  equivalents.  As  a  result  of  this  effort,  the  model,  with  very 
minor  modifications,  has  run  successfully  on  the  following  computers:  CDC-6600, 
CDC-6400,  CYBER-750.  CRAY-1.  Harris-800,  and  FPS-1  array  processor  at¬ 
tached  to  a  VAX  11/780. 

In  addition  to  flexible  resolution,  the  model  can  be  run  with  a  choice  of  rhom- 
boidal  or  triangular  truncation.  For  rhomboidal  truncation,  the  spectral  expan¬ 
sions  are  of  the  form 

M  |m|  +  M  . 

f(A.<£.  a,  t)  =  £  £  f“  (O.  t)  P™  (sin^)  elmA  (18) 

m  =  -M  n  =  |m|  n  n 

with  (2M  +  1)(M  +  1)  spectral  degrees  of  freedom,  and  P™  is  the  associated 
Legendre  function  of  order  m  and  degree  n.  The  transform  grid  consists  of  at 
least  3M  +  1  equally  spaced  longitudes  and  at  least  (5M  +  l)/2  Gaussian  latitudes. 
For  triangular  truncation,  the  expansions  are  of  the  form 

MM  .  . 

f(A.  <t>,  a.  t)  =*  £  £  f™  (a,  t)  P™  (sin#)  eimA  (19) 

m  =  -M  n  »  |mj  n 

2 

with  (M  +  1)  degrees  of  freedom,  and  a  transform  grid  of  at  least  3M  +  1  equally 
spaced  longitudes  and  at  least  (3M  +  l)/2  Gaussian  latitudes. 

Figure  1  shows  the  spectral  domain  for  the  low-  and  high -resolution  versions 
of  the  model.  The  low -resolution  versions  are  rhomboidal  15  (indicated  as  15R) 
and  triangular  21  (2 IT)  with  496  or  484  degrees  of  freedom  respectively.  The  high- 
resolution  versions  are  30R  and  42T  with  1891  and  1849  degrees  of  freedom  re¬ 
spectively. 

For  efficiency  and  convenience,  the  rhomboidal  and  triangular  versions  are 
maintained  as  separate  codes.  The  reason  for  this  is  related  to  the  method  of 
storage  of  the  spectral  coefficients.  For  rhomboidal  truncation,  the  spectral  do¬ 
main  (Figure  1)  can  be  easily  mapped  into  a  square  array  as  shown  in  Figure  2a. 
Each  coefficient  is  referred  to  by  two  indices  (mm,  nj)  which  are  related  to  the 
mathematical  indices  (m,  n)  by 


A 

'  i 


30R 


42  T 


M 


J _ I _ L 

I  2 

nj 

(a) 


12  3  M 


nj 

(b) 


Figure  2.  Arrays  for  Storing  Spectral  Coefficients:  (a)  Rhomboidal 
Truncation,  and  (b)  Triangular  Truncation 


m 


mm  -  1 


(20) 


n  =  nj  +  mm  -  2.  (21) 

The  triangular  truncation  could  also  be  handled  in  this  manner.  However,  it 
would  be  rather  inefficient,  since  almost  half  of  the  array  elements  would  be  zero. 
Instead,  the  coefficients  are  stored  by  columns  in  a  one-dimensional  array 
(Figure  2b).  Each  coefficient  is  now  referred  to  by  a  single  index  ix  obtained 
from  the  relationships 


ix 

=  mm 

+  (nj  *  nj  -  n j)  / 2 

(22) 

m 

=  mm 

-  1 

(23) 

n 

=  nj  - 

1 

(24) 

where  (mm,  nj)  are  the  corresponding  square  array  indices  and  (m,  n)  are  the 
mathematical  indices. 

As  might  be  expected,  an  increase  in  model  resolution  results  in  an  increase 
of  execution  time  and  core  usage.  Table  1  gives  the  computational  requirements 
on  the  CRAY-1  for  various  resolutions  of  the  moist  model.  The  resolution  is 
indicated  as  mRk  or  mTk  where  m  is  the  truncation  wavenumber,  R  or  T  refers  to 
the  type  of  truncation  (rhomboidal  or  triangular),  and  k  is  the  number  of  layers  in 
the  vertical. 

A  comparison  of  15R12  and  21T12  indicates  that,  for  roughly  the  same  number 
of  degrees  of  freedom,  the  triangular  version  runs  somewhat  slower.  Two  reasons 
exist  for  this.  First,  the  2  IT  12  performs  more  computations  because  of  a  larger 
transform  grid  of  64  x  34  (2176  points)  as  compared  to  48  x  40  (1920  points)  for 
18R12.  The  second  reason,  elaborated  upon  in  section  4,  is  that  vectorization  is 
not  as  efficient  for  triangular  truncation. 


i.  \  n<  .ToKizvnoN 

The  first  step  in  transferring  the  model  to  the  C RAY- 1  consisted  of  removing 
all  machine  dependent  structures  (for  example,  multiple  statements  on  a  single 
line)  and  replacing  all  calls  to  machine  language  routines  by  standard  Fortran 
equivalents.  Thus,  with  relatively  little  effort,  a  running  version  of  the  model 
was  put  on  the  CRAY-1.  This  will  be  referred  to  as  "the  original  version  of  the 
model."  From  the  first  three  lines  of  Table  2,  it  can  be  seen  that  the  original 


Table  1.  Computational  Requirements  on  the  CRAY-1  for  ' 
Various  Resolutions  of  the  Moist  Model 


Model  (resolution) 

cpu  time  (aec)/model 
time  step 

core  (K -words) 

15R6 

1. 10 

89 

15R12 

2.02 

137 

2  IT  12 

2.59 

145 

30R12 

9.72 

393 

42T12 

13.68 

408 

version  of  the  model  ran  12-13  times  faster  on  the  CRAY-1  than  on  the  CDC-6600 
with  maximum  compiler  optimization.  However,  the  small  difference  between  the 
scalar  and  vector  execution  times  on  the  CRAY-1  (lines  2  and  3  in  Table  2)  indi¬ 
cates  that  the  original  code  was  not  structured  properly  for  vectorization. 


Table  2.  Effects  of  Restructuring  Code  (Vectorization) 
on  Execution  Time  for  a  24-Hour  Moist  Forecast  Using 
Six  Layers  and  Rhomboidal  M  =  15 


Model  Version  /  hardware  (compiler) 

cpu  time  (sec) 

original  /  CDC-6600 

(opt  «  2) 

1076 

original  /  CRAY-1 

(scalar) 

89 

original  /  CRAY-1 

(vector) 

82 

vectorized  /  CDC-6600 

(opt  =  2) 

679 

vectorized  /  CYBER -750  (opt  =  2) 

225 

vectorized  /  CRA  Y-l 

(scalar) 

49 

vectorized  /  CRAY-1 

(vector) 

27 

To  take  advantage  of  the  vector  processing  capabilities  of  the  CRAY-1,  we 
tried  to  rewrite  the  model  in  a  more  efficient  form.  We  focused  primarily  on 
those  subroutines  that  handle  the  transforms,  since  they  comprise  the  most  time- 
consuming  portion  of  the  code.  The  code  resulting  from  this  restructuring  will  be 
referred  to  as  'the  vectorized  version.  '1  Table  2  shows  that  the  vectorized  version 
of  the  model  runs  more  efficiently  than  the  original  on  both  the  CDC-6600  and  the 
CRAY-1.  The  net  effect  of  restructuring  the  code  (comparing  lines  3  and  7 )  is  a 
speed-up  by  a  factor  of  three.  A  comparison  of  the  original  version  on  the 
CDC-6600  with  the  vectorized  version  on  the  CRA  Y-l  (lines  1  and  7)  shows  an 
overall  speed-up  of  a  factor  of  40. 


Finally,  Table  1  compares  the  relative  efficiencies  of  the  two  different  trunca¬ 
tions.  For  a  comparable  number  of  spectral  degrees  of  freedom  (15R12  compared 
to  21T12),  Tablet  shows  that  the  triangular  truncation  required  28%  more  execu¬ 
tion  time.  As  was  mentioned  in  Section  3,  there  are  two  reasons  for  this,  each 
accounting  for  roughly  half  of  the  additional  time.  First,  the  triangular  version 
has  a  slightly  larger  transform  grid  and  thus  performs  more  grid  point  calcula¬ 
tions.  The  second  reason  is  that  the  triangular  version  is  not  as  amenable  to 
vectorization  because  of  the  structure  of  the  computational  vectors  used  in  the 
Legendre  transform  and  Gaussian  quadrature.  The  difference  is  that  the  rhomb- 
oidal  vectors  are  all  of  a  fixed  length,  while  the  triangular  vectors  are  of  varying 
lengths  and  thus  not  as  amenable  to  efficient  vector  processing.  However,  methods 
exist  to  further  optimize  the  triangular  code  and  reduce  the  gap  between  the  two 
truncations.  Some  of  these  methods  will  be  examined  in  the  near  future. 

5.  REVISED  PREPROCESSING  AND  POSTPROCESSING 

As  in  TR-1  the  1978  FGGE  level  III-A  data  are  still  being  used  as  the  input 
and  verifying  data.  The  terms  "preprocessing"  and  "postprocessing"  were  intro¬ 
duced  in  TR-1  to  describe  the  procedures  required  for  the  transformations  before 
and  after  prognostication,  respectively,  between  the  analysis  (or  data)  grid  and 
the  model  grid. 

Two  revisions  have  been  made  to  the  procedures  described  in  TR-1.  The 
first  is  that  instead  of  the  analyzed  temperature  from  the  FGGE  data,  a  new  set 
of  temperature  values  on  the  mandatory  levels  is  now  used  to  calculate  the  initial 
values  of  temperature  and  geopotential  in  the  model  grid.  These  new  values  are 
obtained  by  solving  a  linear  system  of  equations  that  accord  with  the  least-squares 
principle  and  that  use  the  analyzed  height  data  as  the  input.  The  new  set  of  temp¬ 
erature  values  is  closer  to  being  in  hydrostatic  balance  with  the  set  of  analyzed 
height  values  than  is  the  set  of  analyzed  temperature  values  and  is,  therefore, 
more  suitable  for  input  to  the  model. 

The  second  revision  is  in  the  vertical  interpolation  of  humidity.  Previously, 
values  of  relative  humidity  on  mandatory  pressure  levels  in  the  FGGE  data  set 
were  converted  into  specific  humidity  prior  to  the  vertical  interpolation  step 
where  it  is  assumed  that  the  chosen  humidity  variable  varies  linearly  with  the 
logarithm  of  pressure  between  any  two  consecutive  data  levels.  A  study  by 
Mitchell  and  Yang,  however,  demonstrated  convincingly  that  specific  humidity 

5.  Mitchell,  K. ,  and  Yang,  C.  (1984)  A  Comparison  of  Moisture  Variables  in  the 
Vertical  Interpolation  of  a  4-D  Data  Assimilation  System.  AFGL  Technical 


does  not  vary  linearly  with  the  logarithm  of  pressure.  Thus,  in  conjunction  with 
the  interpolation  and  extrapolation  procedures  currently  used,  relative  humidity 
is  more  suitable  as  the  variable  of  interpolation.  Accordingly,  the  conversion 
from  relative  humidity  to  specific  humidity  in  the  preprocessing  has  been  shifted 
to  follow  the  vertical  interpolation.  Similarly,  in  the  postprocessing,  the  conver¬ 
sion  from  specific  humidity  to  relative  humidity  is  now  done  before  the  vertical 
interpolation. 

In  the  comparison  of  different  resolutions  reported  in  sections  6,  7  and  8,  the 
same  assumptions  and  algorithms  have  been  used  to  perform  both  preprocessing 
and  postprocessing  regardless  of  resolution,  since  the  main  objective  of  the  study 
is  to  assess  the  influence  of  resolution  on  model  forecasts.  In  the  vertical  inter¬ 
polations  between  the  analysis  and  model  grids,  temperature,  wind,  and  relative 
humidity  are  all  assumed  to  vary  linearly  with  the  logarithm  of  pressure  in  the 
layer  bounded  by  any  two  consecutive  levels  where  data  are  available. 

Since  the  model  treats  the  atmosphere  as  a  number  of  discrete  layers,  a 
representative  pressure  level  must  be  defined  for  each  layer.  The  prognostic 
variables  at  these  levels  are  determined  by  interpolating  from  the  analyzed  values 
at  the  adjacent  mandatory  pressure  levels.  It  should  be  noted  that  the  three  phase 
of  the  model  (preprocessing,  forecasting,  and  postprocessing)  do  not  agree  with 
each  other  in  assigning  layer  temperatures  to  pressure  levels. 

A  layer  temperature,  T^,  ^+1,  representing  the  layer  bounded  by  pressures 
Pk  and  pk+i  that  are  at  heights  and  Z^,  respectively,  is  defined  as 


Zk  -  Zk+1 


In  Pk+1 


In  the  preprocessing  phase,  this  temperature  is  assigned  to  a  level  with  pressure 
p'k.  k+»  «‘ven  b* 

ta?k,kH  ‘  »«(l„pk  <26> 


In  the  forecast  phase,  the  temperature  in  the  layer  bounded  by  two  sigma 
interfaces  <7.  ,  <7  is  assigned  to  the  level  with  pressure 


where  p*  is  the  surface  pressure.  Finally,  in  the  postprocessing  phase,  the  same 
layer  temperature  is  assigned  to  the  level  with  pressure 


P  ’ 


k,  k+1 


~  K  ~  K 
rk+l  ak 


*(ln  ffk+l  '  ln  V 


l  /# 


(28) 


Table  3  exhibits  the  numerical  differences  among  these  three  pressures  in 
the  high-resolution  model  where  the  12-layer  structure  of  the  NMC  operational 
model  is  used. 


(t.  >IOI>KI.  PERFORM  VSCK  -  COMPARISON  OF  HIGH-  AND 
1.0  M  -RESOI.I  TION  FORECASTS 

In  the  studies  reported  in  the  following  sections,  "low  resolution"  refers  to 
the  forecast  model  with  rhomboidal  truncation  at  wavenumber  15  in  the  horizontal 
and  6  layers  in  the  vertical  (15R6);  "high  resolution"  refers  to  the  forecast  model 
with  rhomboidal  truncation  at  wavenumber  30  in  the  horizontal  and  12  layers  in  the 
vertical  (30R  12).  The  a  structures  are  summarized  in  Table  4.  A  detailed  com¬ 
parison  of  the  high- resolution  rhomboidal  (30R)  and  triangular  (42T)  will  be  con¬ 
ducted  after  the  triangular  code  is  further  optimized. 

Tables  5  and  f>  present  the  statistics  of  the  24-  and  48-hour  forecast  perfor¬ 
mances  over  periods  of  five  consecutive  days  each  in  January  (15-19  January  1978) 
and  July  (16-20  July  1978).  Each  entry  contains  the  sample  mean  and  the  sample 
standard  deviation  (in  parenthesis)  of  the  global  root-mean-square  (rms)  differences 
(see  TR-1  for  details).  The  standard  deviation  measures  the  representativeness 
of  the  corresponding  mean.  Symbols  A,  F,  S,  and  P  stand  for  analysis,  forecast, 
synthesis,  and  persistence,  respectively.  The  analysis  refers  to  the  fields 
obtained  by  interpolating  FGGE  1 1 1- A  data  from  the  2.5°  -interval  latitude-longi¬ 
tude  grid  to  a  2.  5  longitude  Gaussian  latitude  grid.  These  fields  are  used  both  as 
the  input  and  for  verification.  The  forecast  includes  preprocessing  at  the  initial 
time,  initialization,  prediction,  and  postprocessing  at  the  verifying  time,  while 
the  synthesis  subjects  the  FGGE  data  at  the  verifying  time  to  the  preprocessing 
and  postprocessing  The  difference  between  the  synthesis  and  the  analysis  is 
therefore  a  measure  of  the  direct  impact  of  the  interpolations  and  truncations  in 


Table  3.  Pressure  Values  (mb)  of  Levels  to  Which  Layer  Temperatures 
Are  Assigned  When  the  NMC  12 -Layer  Sigma  Structure  Is  Used.  In  Eqs. 
(26)  and  (28),  the  upper  boundary  is  placed  at  p  ■  10  mb  (p*  ■  1013.25  mb) 


a 

P 

Eq.  (26) 
preprocessing 

Eq.  (27) 
forecasting 

Eq.  (28) 
postprocessing 

0. 

0. 

22.51 

21.02 

23.22 

.050 

50.66 

2 

71.  65 

74.97 

72.06 

.  100 

101.33 

3 
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126. 05 

124.34 

.  150 

151.99 

4 

175.50 

176.89 

175.67 

.200 
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5 

226. 57 

227.65 
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.250 

253.31 

6 

277.49 

278.37 

277.60 

.300 

303.98 

7 

339.85 

341.47 

340.06 

.375 

379.  97 

8 

438.75 

442.22 

439. 18 

.500 

506.63 

9 

877.64 

581.43 

578. 11 

.650 

658.61 

10 

730. 66 

733.67 

731.04 

.800 

810.60 

11 

871. 63 

873.38 

871.85 

.925 

937.26 

12 

974.51 

975.08 

974.58 
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1013.25 


Table  4.  The  a -Structures 
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0. 
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.062240 
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.050 
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.050 
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.  150 

4 

.050 

. 174573 
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.200 

5 

.050 
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.250 

6 

.050 

.274729 

.300 

3  .250 

.369929 

7 

.075 

.337003 

.375 

8 

.  125 

.436433 
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.500 

4  .250 
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9 

.  150 

.573831 
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10 

.  150 
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.800 

.  900 

11 
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.961960 

6  . 100 

. 949686 
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12 

.075 

.962326 

1.000 

1.000 

Table  5.  Performance  (Global  RMS  Difference  and  Standard  Deviations) 
January  1978.  L  =  low  resolution,  and  H  =  high  resolution 


1.34  1.32 

(.11)  (.13) 


6.92 

(.42) 


•Relative  humidity  is  at  the  300-mb  level 


3.43 

(.22) 


101.40 

(11.17) 


RH  (%) 

24h 

48h 

17.96 

21.82 

(1.05) 

(.41) 

16.76 

21.22 

(.45) 

(.53) 

11.94 

11.82 

(.58) 

(.67) 

3.94 

4.02 

(.19) 

(.20) 

10.90 

12.50 

(.64) 

(.  94) 

23.00  26.  S8  | 
(1.02)  (.44)  ! 


25.56  30.78 
(.81)  (.34) 


25.00  29.96 
(.71)  (.63) 


14.  90 

16.00 

(.62) 

(.58) 

3.32 

3.50 

(.30) 

(.14) 

21.34 

(.40) 


24.64 
(.  68) 
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Table  6.  Performance  (Global  RMS  Differences  and  Standard  Deviations)  July  1978. 
L  ■  low  resolution*  and  H  *  high  resolution 


•Relative  humidity  is  at  the  300 -mb  level 


RH  (%) 

24h 

48h 

15.96 

20.  62 

(.63) 

(.31) 

15.46 

20.36 

(.38) 

(.43) 

9.20 

9. 16 

(.24) 

(.31) 

3.44 

3.38 

(.27) 

(.31) 

9.64 

11.90 

(.30) 

(.37) 

89.92 

(5.77) 


22.08  26.  12 
(.58)  (.62) 


19.74  25.40 
(.71)  (.78) 


15.08  15.20 

(.54)  (.39) 


5.34  5.38 

(.23)  (.22) 


21.  10  24.52 
(.69)  (.45) 


26.00  30.06 
(.57)  (.  95) 


26.14  29.78 
(.49)  (.45) 


15.20  15.20 

(.23)  (.23) 


4.84  4.76 

(.30)  (.23) 


20.62  23.74 
(.61)  (.42) 
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the  transformations  between  the  data  and  model  grids.  Persistence  evaluates  the 
difference  of  analyses  between  the  initial  and  verifying  times  (Figure  3). 
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Figure  3.  Schematic  Diagram  of  Model  Evaluation 


The  most  outstanding  difference  between  the  two  resolutions  appears  in  the 
processing  errors  represented  by  the  difference  (S,  A).  The  improved  accuracy 
of  the  high  resolution  is  found  in  all  variables  and  at  all  levels.  The  magnitudes 
of  the  processing  errors  and  the  differences  due  to  resolution  are  similar  in  the 
January  and  July  samples.  When  resolutions  were  varied  both  horizontally  and 
vertically,  it  became  clear  that  most  of  these  improvements  resulted  from  the 


•  ,N 
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horizontal  resolution.  Figure  4,  a  meridional  cross-section  of  the  rms  differences 
between  the  synthesized  and  analyzed  500  mb  heights,  shows  such  an  example.  The 
only  exception  to  this  rule  was  found  in  the  temperature  at  the  250  mb  level  which 
was  affected  more  by  the  vertical  resolution,  as  shown  in  Figure  5. 


LPT I TUOE 


Figure  4.  Meridional  Profile  of  Processing  Error  -  RMS  Error 
of  500  mb  Heights  for  Various  Resolutions  (All  Rhomboidal 
Truncation)  for  FGGE  Data  of  00Z  18  January  1978 


The  main  improvements  in  the  forecasts,  on  the  other  hand,  are  found  in  the 
heights  at  all  three  levels  and  in  the  temperature  at  the  250  mb  level.  In  the  latter 
case,  almost  all  of  the  improvement  in  the  forecasts  appears  to  result  from  the 
reduction  in  the  processing  errors  produced  by  the  increased  vertical  resolution. 
In  the  improved  height  forecasts,  it  is  difficult  to  determine  whether  the  reduction 
in  error  is  due  to  the  reduced  initial  error  or  to  the  reduced  error  growth  rate. 


LATITUDE 


Figure  5.  Meridonial  Profile  of  Processing  Error  -  RMS  Error  of 
250  mb  Temperatures  for  Various  Resolutions  (All  Rhomboidal 
Truncation)  for  FGGE  Data  of  00Z  18  January  1978 


We  will  try  to  answer  this  question  in  the  future  by  running  and  comparing  addi¬ 
tional  resolution  experiments.  One  additional  striking  feature  of  the  forecast 
errors  is  that  the  relative  humidity  forecasts  are  less  skillful  than  persistence. 
Some  of  the  possible  reasons  for  this  will  be  discussed  below. 

Figures  6  and  7  summarize  some  of  the  more  important  aspects  of  the  per¬ 
formances  and  relative  skills  of  the  high-  and  low -resolution  models.  Figure  6 
shows  the  growth  rate  of  the  global  rms  error  of  height  for  the  troposphere 

(average  for  1000  to  200  mb)  for  six  forecasts  (three  in  January  and  three  in  July). 

6 

Persistence  errors  approach  the  climatic  variance  (typical  value  of  about  100  ni) 


6.  Bengstsson,  L.  (1981)  Numerical  prediction  of  atmospheric  blocking:  A  case 
Study.  Tellus  33:19-42. 


within  two  to  three  days.  Beyond  three  days,  the  curve  flattens  out  and  shows 
very  little  increase.  If  we  define  persistence  as  the  lower  limit  of  a  useful  fore¬ 
cast.  then  it  can  be  seen  that  the  low  resolution  (15R6>  height  forecasts  have  no 
skill  beyond  four  days.  The  high-resolution  height  forecasts  (30R12)  are  still 
useful  at  four  days  and  probably  retain  some  skill  at  six  or  seven  days. 

- PERSISTENCE 

- LOW  RES05R6) 

. HIGH  RES  (30RI2) 


0  12  3  4 


FORECAST  DAY 

Figure  6.  Growth  Rate  of  the  Average  1000-200  mb  Height  Errors  for  Six  Fore¬ 
casts  (Three  From  January  1978  and  Three  From  July  1978) 

Figure  7  shows  the  error  growth  rate  for  the  850-mb  relative  humidity  for  the 
same  six  forecasts.  In  contrast  to  the  other  prognostic  variables,  relative 
humidity  forecasts  have  no  skill  at  all  as  measured  by  a  comparison  with  persis¬ 
tence.  Furthermore,  very  little  difference  exists  between  the  errors  of  the  low- 
and  high-resolution  forecasts.  This  indicates  that,  for  reasons  yet  to  be  deter¬ 
mined,  there  is  a  very  strong  tendency  for  the  model -predicted  humidity  to  quickly 
evolve  to  a  certain  preferred  or  "climatological"  state  that  may  differ  substantially 
from  the  Hough -analyzed  humidity.  This  characteristic  seems  to  be  common  to 


all  global  models.  ’  The  other  interesting  point  to  note  in  Figure  7  is  the  very 
small  growth  rate  in  persistence  beyond  day  one.  This  lack  of  variation  in  the 
humidity  analyses  is  probably  due  to  the  coarse  resolution  and  poor  quality  of  the 
Hough  analysis  scheme  as  applied  to  relative  humidity. 
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Figure  7.  Growth  Rate  of  the  Average  850  mb  Relative  Humidity  Errors 
for  Six  Forecasts  (Three  From  January  1978  and  Three  From  July  1978) 


An  example  of  the  low-  and  high-resolution  forecasts  is  shown  in  Figures  8, 

9,  10,  and  11.  Each  figure  includes  four  northern  hemisphere  polar  stereographic 
plots  of  height:  (a),  the  verifying  analysis,  (b)  the  low  resolution  forecast  (15R6), 
(c),  the  high  resolution  forecast  (30R12),  and  (d)  the  difference  between  the  high 
and  low  resolutions.  The  initial  data  consisted  of  the  FGGE  III-A  analysis  for 
00Z  17  January  1978.  Figures  8  and  9  show  the  48-hour  forecasts  of  the  1000  mb 
and  500  mb  heights,  respectively.  At  48  hours,  both  resolutions  can  reasonably 


7.  Sirutis,  J. ,  Miyakoda,  K. ,  and  Ploshay,  J.  (1980)  Moisture  distribution  de¬ 

rived  in  mathematical  models  and  four  dimensional  assimilation,  in 
Atmospheric  Water  Vapor,  Academic  Press,  New  York,  pp.  489-496. 
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Nieminen,  R.  (1983)  Operational  Verification  of  ECMWF  Forecast  Fields  and 
Results  for  1980-1981,  ECMWF  Technical  Report  No.  36. 
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Figure  8a.  48-Hour  Forecast  of  1000  mb  Heights  Beginning  From  00Z  17  January 
1978:  Verifying  Analysis.  Contour  interval  is  50  m 
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Figure  8b.  48-Hour  Forecast  of  1000  mb  Heights  Beginning  From  OOZ  17  January 
1978:  Low  Resolution  (15R6).  Contour  interval  is  50  m 
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Figure  8c.  48-Hour  Forecast  of  1000  mb  Heights  Beginning  From  OOZ  17  January 
1978:  High  Resolution  (30R12).  Contour  interval  is  50  m 


48  HR  FORECAST  HEIGHT  DIFFERENCE  BETWEEN  HIGH 
AND  LOW  RESOLUTIONS 
IOOO  MB 
OOZ  JAN  19, 1978 


Figure  8d.  48-Hour  Forecast  of  1000  mb  Heights  Beginning  From  OOZ  17  January 
1978:  Difference  Between  High  and  Low  Resolutions.  Contour  interval  is  50  m 
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Figure  9a.  48-Hour  Forecast  of  500  mb  Heights  Beginning  From  00Z  17  January 
1978:  Verifying  Analysis.  Contour  interval  is  50  m 
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Figure  9b.  48-Hour  Forecast  of  500  mb  Heights  Beginning  From  00Z  17  January 
1978:  Low  Resolution  (15R6).  Contour  interval  is  50  m 
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Figure  9c.  48-Hour  Forecast  of  500  mb  Heights  Beginning  From  OOZ  17  January 
1978:  High  Resolution  (30R12).  Contour  interval  is  50  m 
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Figure  9d.  48-Hour  Forecast  of  500  mb  Heights  Beginning  From  OOZ  17  January 
1978:  Difference  Between  High  and  Low  Resolutions.  Contour  interval  is  50  m 
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Figure  10a*  96-Hour  Forecast  of  1000  mb  Heights  Beginning  From  00Z  17  January 
1978:  Verifying  Analysis.  Contour  interval  is  50  m 
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Figure  10b.  96-Hour  Forecast  of  1000  mb  Heights  Beginning  From  OOZ  17  January 
1978:  Low  Resolution  (15R6).  Contour  interval  is  50  m 
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Figure  10c.  96-Hour  Forecast  of  1000  mb  Heights  Beginning  From  OOZ  17  January 
1978:  High  Resolution  (30R12).  Contour  interval  is  50  m 
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Figure  lOd.  96-Hour  Forecast  of  1000  mb  Heights  Beginning  From  OOZ  17  January 
1978:  Difference  Between  High  and  Low  Resolutions.  Contour  interval  is  50  m 
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Figure  11a.  96-Hour  Forecast  of  500  mb  Heights  beginning  From  00Z  17  January 
1978:  Verifying  Analysis.  Contour  interval  is  50  m 
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Figure  lib.  96 -Hour  Forecast  of  500  mb  Heights  Beginning  From  OOZ  17  January 
1978:  Low  Resolution  (15R6).  Contour  interval  is  50  m 
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Figure  11c.  96-Hour  Forecast  of  500  mb  Heights  Beginning  From  OOZ  17  January 
1978:  High  Resolution  (30R12).  Contour  interval  is  50  m 
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Figure  lid.  96-Hour  Forecast  of  500  mb  Heights  Beginning  From  OOZ  17  January 
1978:  Difference  Between  High  and  Low  Resolutions.  Contour  interval  is  50  m 


forecast  the  locations  of  the  major  pressure  systems.  In  general,  however,  the 
amplitudes  in  the  low-resolution  forecast  tend  to  be  too  small.  This  is  most 
noticeable  in  the  strong  cyclone  that  has  developed  over  the  mid-Pacific  (at  both 
1000  and  500  mb).  Also,  the  low -resolution  forecast  cannot  adequately  capture 
smaller  scale  features  such  as  fronts. 

Figures  10  and  11  show  the  96-hour  forecasts  of  the  1000  mb  and  500  mb  heights, 
respectively.  Here,  the  differences  in  resolution  become  much  more  apparent. 

The  low -resolution  forecast  has  lost  most  of  its  predictive  skill  (also  see  Figure 
12)and  begins  to  suffer  from  both  amplitude  and  phase  errors.  This  is  especially 
noticeable  in  the  500  mb  forecast  and  in  the  differences  between  the  high  and  low 
resolution  forecasts  (Figure  lid). 

Finally,  in  Figure  12,  we  show  the  growth  rates  of  the  height  errors  at 
(a)  1000  mb,  (b)  500  mb,  (c)  250  mb,  and  (d)  average  1000-200  mb.  In  general, 
the  difference  in  skill  between  the  low-  and  high-resolution  forecasts  increases 
with  time  and  altitude.  The  extremely  poor  performance  of  the  low-resolution 
forecast  at  250  mb  is  caused  primarily  by  the  lack  of  adequate  vertical  resolution 
in  the  upper  troposphere  and  lower  stratosphere. 


7.  REPRESENTATION  OF  TOPOGRAPHY 

The  basic  source  of  information  on  the  surface  topography  is  a  set  of  height 
values  on  the  2.5°  -interval  latitude-longitude  coordinates  furnished  by  NMC  as  a 
part  of  the  fixed  field  data  in  the  FGGE  III-A  data  set.  It  is  assumed  that  the 
value  at  each  grid  point  has  been  obtained  through  an  averaging  process  and  repre¬ 
sents  the  mean  terrain  height  of  the  area  over  which  the  average  has  been  defined. 
The  network  of  grid  points  on  which  various  calculations  are  performed  in  a 
spectral  model  is  generally  different,  from  the  network  of  grid  points  in  the  basic 
data  set.  Consequently,  a  question  arises  about  how  the  representation  of  the 
surface  topography  in  a  spectral  model  of  specified  resolution  should  be  generated 
from  the  basic  data. 

The  theory  of  linear  transformations  shows  that  any  set  of  real-valued  data 

(X  (A..,  6t  ),  j  =  1,  ...  ,  J;  1  1  1,  ...  L>)  such  as  topography  may  be  represented 
J 

through  the  Fourier  transform  by  an  analytic  function 
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Growth  Rate  of  Height  Errors  for  Forecasts  From  00Z  17  January  1978:  (a)  1000  mb,  (b)  500  mb. 
and  (d)  Average  1000-200  mb 


such  that  X  (Xj,  <bt )  =  X  (Xj,  #/ )  at  all  j  and  l .  Here,  M  =  J/2  (when  J  is  even) 
or  (J- 1 '/ 2  (when  J  is  odd),  and  Qm  (sin#)  is  a  polynomial  of  degree  (L  - 1)  in  sin# 
for  each  m  such  that  Q_m  and  Qm  form  a  complex  conjugate  pair. 

The  analytic  function  X(X,  #),  on  the  other  hand,  may  be  represented  in 
terms  of  spherical  harmonics. 


X(X.  #)  = 


E  E 

m  =  -M  n=  |m| 


P™  (sin#)  el 


in  which 


x"  -  /  Qm<y)  p“  d* 


where 


/,  p>’  p™  w 


for  all  m,  y  -  sin#,  and  inf),  is  the  Kronecker  delta.  The  integral  in  Eq.  (31) 
will  be  referred  to  as  "the  Legendre  transform."  The  integrand  in  Eq.  (31)  is 
seen  to  be 


Qm(y)  (y) 


a  polynomial  in  y  of  degree  (L  +  n  -  1 ) 

when  m  is  even  ; 

2  1/2 

(1  -  y  )  times  a  polynomial  in  y  of 
degree  (L  +  n  -  2 )  when  m  is  odd. 


In  practice,  the  wave  amplitudes  at  a  given  latitude  4>i  are  first  obtained  from 
(X  (Xj,  #^  ),  j  =1,  ....  J)  using  the  Fourier  transform 


*,)  ■  E  X<A,.  *L  > 
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and  then  the  integral  in  Eq.  (31)  is  replaced  by  a  quadrature  which  may  be  written 


Lj 

X^  =  £  Q  (y  .  )  P*” (y ,  )  w. 

n  m  •'l  n  *  t  * 
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where  =  sin^  and  w^  is  the  weight  associated  with  .  The  reconstructed 
topography  becomes 


X  (A,  4>)  = 


m=-M 


P“  (y)  eim* 


in  which  MQ  and  NQ  define  the  range  of  truncation.  The  error  of  spectral  synthe¬ 
sis,  X(Aj,  <p^)  -  X(Aj,  4> £  ),  is  seen  to  arise  from  two  possible  sources,  one 
in  the  quadrature  [Eq.  (33)]  and  the  other  in  the  truncation  [Eq.  (34)]. 

The  impacts  of  quadrature  and  truncation  on  the  representation  of  topography 
were  examined  by  comparing  various  truncations  in  Eq.  (34)  and  by  using  the 
trapezoidal  and  Gauss-Legendre  (G-L)  quadratures  in  Eq.  (33).  The  trapezoidal 
quadrature  assumes  that  the  variation  of  the  integrand  between  any  two  consecu¬ 
tive  data  points  is  linear,  while  the  G-L  quadrature  with  N  Gaussian  coordinates 
evaluates  the  integral  accurate  to  the  polynomial  of  degree  2N  -  1.  The  former 
adapts  itself  to  given  data  distributions,  while  the  latter  requires  data  on  specific 
coordinates  dictated  by  the  degree  of  accuracy.  The  impact  was  measured  by  the 
global  rms  error,  E,  defined  by 


££{x  (Ay  <pt)  -  X (Aj,  *4)}  2Wj  1 

i  j  n 


where  w.  =  w,  in  view  of  the  fact  that  grid  points  are  uniformly  spaced 
around  a  latitude  circle.  When  using  the  G-L  quadrature,  if  the  original  terrain 
data  are  not  on  Gaussian  latitudes,  it  is  necessary  to  estimate  the  values  of  terrain 
height  on  the  Guassian  latitude  at  the  corresponding  longitudes.  This  was  done  by 
a  linear  interpolation  in  sin^  using  the  two  consecutive  data  latitudes  that  sur¬ 
round  the  Gaussian  latitude  of  concern.  The  procedures  encountered  and  quanti¬ 
ties  defined  in  a  full  cycle  of  transformation  are  illustrated  and  designated  in 
Figures  13a,  13b,  and  13c. 

In  view  of  the  limited  invertibility  of  the  numerical  quadratures  for  the 
Legendre  transform,  three  measures  of  differences  were  used  to  characterize 
various  aspects  of  the  impacts  of  quadrature  and  truncation  on  the  process  of 
transformation.  Ejis  the  conventional  error  of  synthesis  and  measures  the  com¬ 
bined  impact  of  quadrature  and  truncation.  Eg,  on  the  other  hand,  measures  the 
error  incurred  in  completing  a  full  circle  of  the  transformation  at  a  fixed  trunca¬ 
tion  and  represents  the  error  entirely  due  to  quadrature.  Both  E^  and  Eg  are 
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Figure  13b.  Definitions  of  Quantities,  Processes,  and  Errors  in  Representation  of  Topography  Using 
Gaussian  Quadrature,  Data  on  Equally  Spaced  Latitudes 


jaussian 


defined  in  the  physical  domain.  The  third  measure.  Eg,  is  an  equivalent  of  Eg 
in  the  spectral  domain  and  will  be  identical  in  the  total  magnitude  with  Eg  when 
the  spectral  transform  is  exactly  invertible.  The  composition  of  Eg  reveals  the 
spectral  distribution  of  the  error  of  transformation  due  to  quadrature.  Addition¬ 
ally,  in  order  to  separate  the  effect  of  the  required  interpolation  from  the  original 
data  to  the  Gaussian  latitudes  from  the  effect  of  the  G-L  quadrature,  similar 
measures  of  differences  are  defined  with  reference  to  the  estimates  on  the 
Gaussian  latitudes  and  are  denoted  by  primes  (Figure  13c). 

These  processes  of  transformation  were  applied  to  two  fields  of  topography. 
The  so-called  unsmoothed  terrain  is  the  one  given  in  the  FGGE  III-A  data  set. 

The  so-called  smoothed  terrain  has  been  obtained  by  subjecting  the  unsmoothed 
terrain  to  a  nine-point  smoother  twice.  The  smoother  is  a  product  of  two  three- 
point  smoothers,  one  along  the  zonal  direction  and  the  other  in  the  meridonial 
direction.  It  may  be  represented  by  a  linear  operator  W:  X  =  W  j  X^  defined 

by  W  (i,  k:  j,  i  )  =  w(  i :  j)  w  (k  :  i  ) 
where 


w  (i  :  j  )  = 


1/2  if  i  =  j 
1/4  if  i  =  j  +_  1 
0  otherwise 


and  similarly  for  w(k  :  i  ). 

Such  a  smoother  has  a  progressively  strong  damping  effect  toward  short  waves 
and  completely  eliminates  the  two-grid  interval  wave.  This  is  clearly  evident 
in  Table  7  which  presents  the  amounts  of  variance  contributed  by  various  spectral 
ranges  in  both  the  unsmoothed  and  smoothed  terrains.  These  spectra  were 
obtained  using  the  G-L  quadrature.  The  meridional  profiles  of  three  parameters 
of  these  terrain  fields  are  shown  in  Figures  14a  and  14b.  Here,  DM,  represented 


Table  7.  Amounts  of  Variance  in  Various  Spectral  Ranges  in  the  Unsmoothed 
and  Smoothed  Terrains  (Units  of  m") 
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Figure  14a.  Meridional  Profiles  of  Terrain  Parameters:  DM-Zonal 
Mean.  DE-Zonal  Mean  Plus  Standard  Deviation,  DMAX-Maximum  Value 
for  Unsmoothed  Terrain 


by  a  solid  curve,  is  the  zonal  mean;  DE,  represented  by  a  dotted  curve,  is  the 
sum  of  the  zonal  mean  and  the  standard  deviation  of  the  field  on  the  latitude ;  and 
DMAX,  represented  by  a  solid  curve  with  crosses,  is  the  maximum  value  of  the 
field  on  the  latitude. 

Tables  8 A  and  8B  tabulate  the  global  rms  of  both  the  error  in  synthesis  and 
the  error  in  reproduction  for  both  terrain  fields  at  various  truncations  using 
either  the  trapezoidal  or  G-L  quadratures  with  the  76  Gaussian  latitudes.  The  76 
Gaussian  latitudes  form  the  smallest  set  required  for  the  spectral  model  with  the 
rhomboidal  30  truncation.  The  difference  in  characteristics  between  the  error 
in  synthesis,  Ej,  and  the  error  in  reproduction.  Eg,  is  readily  seen  in  the 
opposite  trends  of  the  variation  of  magnitude  with  truncation  range.  E  decreases 
with  widening  of  truncation  range  as  more  of  the  spectral  components  in  the 
original  fields  are  included.  Eg,  on  the  other  hand,  increases  with  widening  of 
truncation  range  because  of  an  increase  in  the  number  of  polynomials  whose 
degrees  exceed  the  highest  degree  resolvable  by  the  truncation. 

Both  tables  support  the  preference  of  the  G-L  over  the  trapezoidal  quadra¬ 
ture.  Although  the  trapezoidal  quadrature  produces  a  slightly  smaller  rms  error 
in  synthesis,  Ej,  than  the  G-L  quadrature  for  the  unsmoothed  terrain  up  to  the 
rhomboidal  50  truncation,  this  small  edge  is  more  than  compensated  by  disadvan¬ 
tages  found  in  other  aspects.  For  the  smoothed  terrain,  the  G-L  quadrature 
produces  a  smaller  Ejfor  all  truncations.  More  significantly,  the  rms  error  in 
reproduction  increases  more  rapidly  with  the  trapezoidal  quadrature  beyond  the 
rhomboidal  40  truncation  in  the  unsmoothed  terrain  and  across  the  entire  range 
in  the  smoothed  terrain.  It  should  also  be  noted  that  the  global  mean  error  is 
negligible  in  all  cases. 

Further  support  for  favoring  the  G-L  quadrature  is  provided  by  Table 
9  which  presents  Ej  and  Eg  (see  Figure  13c  for  definition)  for  the  two  terrain 
fields.  This  table  compiles  the  errors  of  transformation  if  the  terrain  fields  had 
been  available  or  defined  on  the  Gaussian  latitudes.  The  differences  between  the 
corresponding  quantities  in  Tables  8a,  8b,  and  9  represent  the  effect  due  to  the 
extra  step  of  interpolation  from  the  2.  S’  -interval  latitudes  to  the  Gaussian  lati¬ 
tudes  required  in  obtaining  the  measures  Ej  and  Eg  using  the  G-L  quadrature.  In 
terms  of  rms  error  in  synthesis,  this  amounts  to  approximately  30  ~  36  m  in  the 
unsmoothed  terrain  and  3  ~  5  m  in  the  smoothed  terrain.  The  absence  of  the  inter¬ 
polation  step  in  the  calculation  of  Eg  brings  forth  the  complete  invertibility  of  the 
G-L  quadrature  as  long  as  the  truncation  range  does  not  exceed  that  specified  by 
the  number  of  Gaussian  latitudes  employed.  The  76  Gaussian  latitudes  should 
reproduce  exactly  up  to  the  rhomboidal  37  truncation,  beyond  which  the  error  in 
reproduction  should  increase  with  further  widening  of  truncation  range.  Table  9 
bears  witness  to  these  theoretical  inferences.  In  fact,  the  values  of  Eg  at 
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rhomboidal  60  and  70  truncations  exceed  those  of  Eg. 


Table  8a.  Root -Mean -Squares  of  the  Error  in  Synthesis  (Ej)  and  Error  in  Repro¬ 
duction  (E2)  With  Different  Rhomboidal  Truncations  of  Unsmoothed  Terrain  (Units 
of  m) 


Error 

Quadrature 

Rhomboidal  Truncation 

15 

24 

— 

30 

40 

50 

60 

70 

E1 

T 

272.44 

206.32 

175.08 

137.79 

103.20 

88.66 

96.49 

G-L 

272.48 

206.70 

176. 12 

141. 16 

109.92 

88.46 

78.09 

E2 

T 

7. 16 

12.95 

18.21 

25.08 

54.88 

112.39 

174.76 

G-L 

6.96 

14.64 

19.92 

26. 19 

32.03 

38.96 

50.88 

Table  8b.  Root-Mean-Squares  of  the  Error  in  Synthesis  (E^)  and  Error  in  Repro¬ 
duction  (E2)  With  Different  Rhomboidal  Truncations  of  Smoothed  Terrain  (Units  ofm) 


1  - - - - - - - -  1 

Error 

Quadrature 

Rhomboidal  Truncation 

15 

24 

30 

40 

50 

60 

70 

E1 

T 

123.19 

55.61 

37.58 

29.72 

30.24 

33.81 

37.74 

G-L 

123. 14 

55. 14 

35.88 

24.46 

19.46 

16.58 

14.81 

E2 

T 

7. 10 

12.54 

17.05 

24.83 

35.06 

47.17 

57.39 

G-L 

5.48 

8.24 

8.91 

9.44 

9.63 

9.65 

9.65 

Comparisons  of  the  rms  of  errors  in  Tables  8a,  8b,  and  9  of  the  unsmoothed 


and  smoothed  terrain  fields  show  the  large  contribution  made  by  the  smoothing 

,  / 


operation  in  reducing  both  types  of  error,  except  in  Eg  at  truncation  ranges 
where  complete  invertibility  exists. 


Tables  10A  and  10B  summarize  the  statistics  of  the  errors  in  reproduction 


in  the  spectral  domain.  Eg  and  Eg  for  both  terrain  fields.  Comparison  with  the 
corresponding  quantities  in  the  physical  domain,  shows  the  G-L  quadrature  pro¬ 
duces  smaller  differences  between  the  two  domains  than  does  the  trapezoidal 
quadrature,  while  both  quadratures  exhibit  similar  characteristics  in  the  varia¬ 
tions  of  magnitude  with  truncation  range  as  observed  in  Eg  and  Eg. 

The  spectral  power  of  Eg  is  found  in  seven  subranges  in  each  of  the  zonal  and 
meridional  spectra.  The  subranges  are  defined  in  terms  of  an  integer  K,  which 
is  one  plus  the  zonal  wave  number,  m,  for  the  zonal  spectrum  and  one  plus  the 
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Table  9.  Root-Mean-Squares  of  the  Error  in  Synthesis  (Ei)  and  Error  in  Repro¬ 
duction  (E2)  With  Different  Rhomboidal  Truncations  With  Reference  to  the  Given 
Values  on  the  Gaussian  Latitudes  of  the  Unsmoothed  and  Smoothed  Terrain  (Units 
of  m). 


Rhomboidal  Truncation 


■ 

15 

24 

30 

40 

50 

60 

70 

unsmoothed 

238.98 

170.96 

140. 18 

106.16 

77.36 

58.34 

50.21 

smoothed 

117.61 

51.25 

32.33 

20.80 

15.42 

12.01 

9.  68 

unsmoothed 

R2 

R2 

R2 

-5 

9  x  10 

15.65 

43.04 

81.06 

smoothed 

R2 

R2 

R2 

2  x  10"6 

0.19 

0.83 

1.79 

Rg  *  2  x  10'7  is  considered  to  be  round-off  error 

difference,  n  -  m,  between  the  degree  and  order  of  a  spherical  harmonic  compo¬ 
nent  P™  for  the  meridional  spectrum.  The  seven  subranges  are  given  in  Table  11 
[  note  that  these  subranges  can  be  identified  as  horizontal  or  vertical  lines, 
respectively,  in  Figure  2a  as  stated  in  Eqs.  (20)  and  (21)].  Tables  12A  and  12B 

summarize  the  power  distributions  for  the  trapezoidal  quadrature  and  the  G-L 

2 

quadrature,  respectively.  The  total  powers  are  given  in  units  of  m  and  the 
fractional  powers  in  percent. 

First,  with  the  use  of  the  trapezoidal  quadrature,  the  spectral  distributions 
of  the  error  divide  the  truncation  ranges  considered  into  two  groups,  those  with 
rhomboidal  truncation  lower  than  or  higher  than  40.  All  truncations  up  to 
rhomboidal  40  show  the  same  spectral  distributions  in  both  the  zonal  and  meri¬ 
dional  spectra.  The  magnitudes  of  total  power  are  approximately  equal  for  the 
unsmoothed  and  smoothed  terrains.  All  of  the  error  in  the  zonal  spectrum  is 
contained  in  the  zonal  mean,  while  the  error  in  the  meridional  spectrum  increases 
towards  the  short-wave  end.  On  the  other  hand,  for  truncations  exceeding  rhom¬ 
boidal  40,  a  large  error  is  found  at  the  short-wave  end  of  the  zonal  spectrum  for 
the  unsmoothed  terrain.  This  error  arises  from  the  contributions  of  the  short 
waves  in  the  unsmoothed  terrain  that  are  effectively  damped  out  by  the  smoothing 
operator.  This  is  evident  from  the  lack  of  any  corresponding  short-wave  error  in 
the  smoothed  terrain  for  which  the  error  is  all  in  the  zonal  mean.  Also,  the  total 
power  of  the  smoothed  terrain  is  about  the  same  as  power  of  the  zonal  mean  of  the 
unsmoothed  terrain.  The  distributions  of  the  meridional  spectra  of  the  two  terrain 
fields  are  similar,  although  the  unsmoothed  spectrum  is  more  skewed  towards  the 
short-wave  end. 


Table  10a.  The  Square  Roots  of  the  Power  of  the  Error  in  Reproduction  in  the 
Spectral  Domain  (E3  or  Eg)  in  the  Unsmoothed  Terrain  (Unit  ■  m) 


Quadrature 

Error 

Rhomboidal  T  runcation 

15 

24 

30 

40 

50 

60 

70 

T 

E3 

7.02 

12.36 

16.99 

22.40 

42.50 

82.09 

125.60 

G-L 

E3 

6.96 

14.64 

19.  92 

26.20 

31.72 

37.09 

49. 10 

R3 

R3 

R3 

11.30 

30.77 

57.61 

_7 

Rj  *  2  x  10  is  considered  to  be  round-off  error 


Table  10b.  The  Square  Roots  of  the  Power  of  the  Error  in  Reproduction  in  the 
Spectral  Domain  (E3  or  E3)  in  the  Smoothed  Terrain  (Unit  ■  m) 


Quadrature 

Error 

Rhomboidal  Truncation 

15 

24 

30 

40 

50 

60 

70 

T 

E3 

6.96 

11.97 

15.91 

22. 18 

29.87 

38.45 

45.64 

m 

E3 

5.48 

8.24 

8.91 

9.44 

9.63 

9.68 

9.75 

UU 

R3 

R3 

R3 

mm 

0. 19 

0.87 

1.74 

R3  s  2  x  10'7  is  considered  to  be  round-off  error 


Unlike  the  trapezoidal  quadrature,  the  G-L  quadrature  exhibits  considerable 
differences  in  both  the  total  magnitudes  and  the  spectral  distributions  of  the  error 
for  the  two  terrains  for  all  truncation  ranges.  For  the  unsmoothed  terrain  with 
truncations  less  than  40,  the  G-L  quadrature  produces  errors  that  are  larger 
than  the  trapezoidal  quadrature.  For  truncations  beyond  rhomboidal  40,  the  G-L 
errors  are  significantly  less  than  the  trapezoidal  errors.  For  the  smoothed 
terrain,  however,  the  G-L  quadrature  produces  errors  that  are  much  less  than 
those  of  the  trapezoidal  quadrature.  These  errors  also  increase  much  more 
slowly  with  the  widening  of  truncation  range  than  all  others.  The  spectral  distri¬ 
butions  tend  to  be  mound-shaped  for  the  smoothed  terrain  for  both  the  zonal  and 
meridional  spectra.  For  the  unsmoothed  terrain,  the  distributions  are  skewed 
towards  the  short-wave  end. 

The  influence  of  the  differences  in  the  terrain  fields  on  model  performance 
was  assessed  in  terms  of  the  global  rms  errors  of  the  height  forecasts  of  the 
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Table  12b.  Total  Power  (in  m  )  and  Spectral  Composition  (in  Percent)  of  Eg  for 
Gauss- Legendre  Quadrature  (L  Zonal,  M=  Meridional) 


Rhomboldal  Truncation 


15 

24 

30 

40 

50 

60 

70 

a 

i 

U 

S 

U 

S 

U 

S 

U 

S 

U 

s 

U 

1 

Total 

Power 

46 

m 

214 

68 

397 

79 

686 

89 

1011 

1376 

94 

2411 

m 

1 

L 

8 

n 

4 

71 

3 

8 

3 

7 

3 

8 

3 

8 

2 

7 

M 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

L 

38 

46 

24 

m 

17 

a 

16 

38 

■9 

39 

13 

39 

8 

38 

2 

M 

6 

8 

2 

H 

1 

1 

1 

1 

■ 

3 

0 

3 

0 

3 

I, 

26 

mm 

23 

28 

16 

27 

KOI 

■>‘Si 

ES 

28 

16 

28 

10 

28 

3 

M 

27 

I 

9 

18 

G 

16 

i 

3 

14 

3 

14 

6 

13 

L 

28 

18 

19 

14 

16 

14 

15 

15 

14 

14 

12 

14 

8 

15 

4 

M 

67 

GO 

23 

29 

15 

26 

10 

23 

8 

23 

7 

22 

11 

22 

L 

0 

0 

31 

29 

0 

23 

9 

19 

8 

16 

12 

9 

r> 

M 

0 

0 

66 

44 

39 

30 

35 

22 

34 

19 

1 

16 

33 

L 

D 

0 

0 

0 

19 

2 

15 

2 

13 

11 

2 

9 

2 

6 

M 

H 

0 

0 

0 

33 

15 

22 

12 

16 

11 

14 

12 

11 

11 

L 

0 

0 

0 

0 

0 

12 

0 

19 

1 

30 

0 

52 

1 

7 

M 

H 

0 

0 

0 

0 

0 

33 

1! 

50 

15 

ra 

16 

56 

17 

mandatory  pressure  levels.  For  this  purpose,  a  comparison  was  made  among 
six  72-hour  forecasts,  three  forecasts  with  each  terrain  field  (unsmoothed  and 
smoothed)  beginning  from  00Z  on  15,  16,  and  17  January  1978.  The  results  are 
summarized  in  Figures  15  and  16. 

Figure  15  shows  the  vertical  profile  of  the  group  means  of  the  global  rms  of 
(1)  the  processing  error  at  the  initial  time,  and  (2)  the  forecast  errors  at  days  one, 
two,  and  three.  It  is  quite  obvious  that  little  difference  exists  in  both  the  magni¬ 
tudes  and  shapes  of  the  profiles  as  a  result  of  differences  in  the  terrain  fields.  It 
is  also  clear  that  the  processing  error  constitutes  a  small  fraction  of  the  forecast 
errors.  No  discernible  difference  exists  in  the  forecast  errors  that  could  be 
ascribed  to  the  difference  in  processing  errors  at  the  initial  time. 

To  probe  further  into  the  relationship  between  the  initial  processing  errors 
and  the  forecast  errors,  the  group  means  and  standard  deviations  of  the 
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differences  between  the  global  rms  errors  for  forecasts  with  the  two  terrain 
fields  were  calculated.  The  results  are  shown  in  Figure  16,  where  a  dot  repre¬ 
sents  the  group  mean  and  the  width  of  the  line  across  the  dot  represents  twice  the 
group  standard  deviation.  Here,  a  positive  value  indicates  a  smaller  error  with 
the  smoothed  topography  and  vice  versa.  From  Figure  16A,  it  is  apparent  that 
use  of  the  smoothed  terrain  reduced  the  processing  error  at  all  levels  except  the 
top  two  (50  mb  and  70  mb).  However,  no  significant  trace  of  this  improvement 
appeared  in  the  forecasts.  The  differences  in  the  forecast  errors  (Figures  16b, 

16c.  and  16d)  were  smaller  than  the  differences  in  the  processing  errors  and  were 
much  smaller  than  the  forecast  errors  themselves.  Furthermore,  these  differ¬ 
ences  were  not  statistically  significant. 

On  the  basis  of  these  findings,  it  can  be  concluded  that,  in  using  a  spectral 
model  for  simulating  and  predicting  the  global  circulation;  (1)  the  terrain  is  best 
defined  on  the  Gaussian  latitudes  of  the  forecast  model;  (2)  the  Gauss- Legendre 
quadrature  is  better  than  the  trapezoidal  quadrature  in  the  computation  of  the 
transforms;  and  (3)  the  smoothed  terrain  is  preferable  as  the  model  terrain. 

Thus,  the  terrain  to  be  used  in  the  high- resolution  model  is  defined  as  the  set  of 
the  spherical  harmonic  coefficients  obtained  from  the  original  FGGE  data  by  first 
passing  them  through  the  9-point  smoother  twice,  interpolating  the  results 
linearly  in  sin0  to  the  76  Gaussian  latitudes,  and  then  transforming  them  into 
spectral  coefficients  at  the  rhomboidal  30  truncation.  The  model  terrain  is 
thereby  uniquely  defined  in  both  the  physical  and  spectral  domains.  For  consis¬ 
tency,  this  terrain  field  is  also  used  in  the  preprocessing  and  postprocessing 
steps. 

8.  EFFECTS  OF  INITIALIZATION 

"Initialization"  is  the  term  generally  used  for  the  procedure  through 
which  unwanted  meteorological  noise  is  suppressed.  This  noise  usually  appears 
as  spurious  gravity  waves  that,  if  left  unchecked,  could  contaminate  and  possibly 
destroy  the  forecast.  These  waves  are  generated  by  certain  imbalances  between 
the  mass  and  motion  fields.  There  are  three  possible  sources  for  these  imbalan¬ 
ces  :  (1)  errors  in  the  observations,  (2)  errors  in  the  analysis  scheme  interpola¬ 
ting  observations  to  grid  points,  and  (3)  the  forecast  model  itself,  since  it  is  a 
discrete  representation  of  a  continuous  fluid.  While  there  is  general  agreement 
within  the  numerical  weather  prediction  community  on  the  need  for  initialization, 
it  is  not  clear  what  form  initialization  should  take  or  to  what  extent  it  is  necessary. 

To  assess  the  impact  of  initialization  on  the  model,  several  high-resolution 
(30R12)  forecasts  using  different  initialization  techniques  have  been  run  and  com- 


pared.  The  methods  that  were  tested  fall  into  two  broad  categories  :  (1)  those  in 
which  the  initial  data  are  filtered  and  balanced  (for  example,  normal  mode  methods) 
and  (2)  those  where  the  spurious  waves  are  controlled  during  the  forecast  by  the 
numerics  and/or  the  parameterized  physics  of  the  model.  Specifically,  three 
different  methods  were  tested:  (1)  no  initialization  at  all,  (2)  linear  damping  of 
divergence  as  described  by  Bourke  and  referred  to  as  "divergence  dissipation," 
and  (3)  the  MachenhauerlO  nonlinear  normal  mode  initialization  process.  All 
three  methods  were  tested  on  forecasts  run  from  data  produced  by  two  different 
analysis  schemes,  the  1978  FGGE  III-A  Hough  analysis (denoted  as  HO)  and  the 
1979  FGGE  III-A  optimal  interpolation  (OI)  analysis.  ^  It  should  be  noted  that  the 
Hough  analysis  implicitly  includes  linear  normal  mode  initialization  since  the 
amplitudes  of  the  gravity  modes  are  set  to  zero  by  the  analysis  scheme.  In  view 
of  this,  it  was  felt  that  Hough  based  forecasts  alone  were  not  sufficient  to  properly 
assess  the  effects  of  initialization,  and  thus  OI  based  forecasts  were  also  run. 

The  forecasts  that  were  run  for  the  six  different  combinations  of  analysis  and 
initialization  and  their  acronyms  are  summarized  in  Table  13. 


Table  13.  Summary  and  Acronyms  of  Different  Combinations 
of  Analysis  and  Initialization  That  Were  Tested 


Initialization 

Method 

Analysis  Scheme 

Hough  (HO) 

Optimal  Interpolation  (OI) 

No  initialization 
(NO) 

HO-NO 

OI-NO 

Divergence  dissipation 
(DD) 

HO-DD 

OI-DD 

Nonlinear  normal  mode 
(NM) 

HO-NM 

OI-NM 

9.  Bourke,  W.  (1974)  A  multi-level  spectral  model.  I.  Formulation  and  hemi¬ 

spheric  integrations,  Mon.  Wea,  Rev,  102:687-701. 

10.  Machenhauer,  B.  (1977)  On  the  dynamics  of  gravity  oscillations  in  a  shallow 

water  model  with  applications  to  normal  mode  initialization,  Beit.  Phys. 
Atmos.  50:253-271 

11.  Flattery,  T.  (1971)  Spectral  models  for  global  analysis  and  forecasting, 

Proc.  Sixth  AWS  Tech.  Exchange  Conf. .  U.S.  Naval  Academy,  Annapolis, 
Md. ,  21-24  September  1970,  AWS-TR-242:42-54. 

12.  McPherson,  R.  D. ,  Bergman,  K.H.,  Kistler,  R.E.,  Rasch,  G.E.,  and 

Gordon,  D.S.  (1979)  The  NMC  operational  global  data  assimilation  system, 
Mon.  Wea.  Rev.  107:1445-1461. 


In  the  no-initialization  cases  (NO),  no  initialization  procedure  is  used 

(except  implicitly  in  the  Hough  analysis  as  noted  above).  The  gravity  waves  are 

controlled  in  the  forecast  model  by  the  semi-imp  licit  time  scheme  that  is  known 

1 3 

to  damp  gravity  waves.  A  certain  amount  of  dissipation  is  also  provided  by  the 

4 

parameterized  subgrid  scale  diffusion,  which  is  a  linear  horizontal  V  operator 
applied  to  vorticity,  divergence,  temperature,  and  specific  humidity.  The  values 
of  the  diffusion  coefficients  are  5  x  10*®  m4  s  *  for  divergence  (applied  to  all 
spectral  components)  and  1  x  10*®  m*  s  *  for  the  other  variables  (applied  only  to 
the  upper  half  of  the  rhomboid,  n  >  M). 

For  divergence  d.ssipation  (DD),  the  initial  values  of  all  spectral  components 
of  divergence  are  set  to  zero.  The  subgrid  scale  diffusion  of  divergence  is 
replaced  by  a  linear  damping  term  of  the  form  -KD  where  the  damping  coefficient, 
K,  is  given  by 

5  x  10'4  exp  (-0.  19188  t)  s'*  t  <  24  h  (3?) 

5  x  10'®  s'*  t  >  24  h 

where  t  is  given  in  hours. 

Diffusion  is  applied  to  vorticity,  temperature,  and  specific  humidity  as  described 
above. 

For  the  nonlinear  normal  mode  initialization  cases  (NM),  the  analyzed  data 
are  filtered  and  balanced  according  to  Machenhauer's10  scheme.  In  this  proce¬ 
dure,  the  data  are  projected  onto  the  normal  modes  of  the  model  and  separated 
into  the  rotational  or  Rossby  modes  and  the  gravity  modes.  The  gravity  modes 
are  adjusted  in  such  a  way  that  the  linear  tendency  terms  are  balanced  by  the 
adiabatic  nonlinear  advective  terms,  thus  resulting  in  a  net  initial  gravity  wave 
tendency  of  zero.  The  procedure  consists  of  two  iterations  applied  only  to  the 

14 

first  four  vertical  modes.  Further  details  of  the  scheme  are  given  by  Ballish. 

In  the  NM  forecasts,  subgrid  scale  diffusion  is  included  as  in  the  NO  cases. 

To  assess  the  impact  of  initialization  on  the  forecasts,  the  discussion  will 
focus  on  the  changes  in  time  of  the  divergence  field  because  it  is  the  quantity  that 
is  most  dramatically  affected  by  initialization.  In  terms  of  the  effect  of  non¬ 
linear  normal  mode  initialisation  on  the  analysed  data,  the  results  were  as 
expected:  very  little  change  in  the  Hough  data  but  substantial  change  in  the  Ol 


13.  Kurihara,  Y.  (1965)  On  the  use  of  implicit  and  iterative  methods  for  the  time 

integration  of  the  wave  equation,  Mon.  Wea.  Rev.  93:33-46. 

14.  Ballish,  A.B.  (1980)  Initialization  Theory  and  Application  to  the  NMC  Spectral 

Model,  Ph.  D.  thesis,  U.  of  Maryland. 
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Figure  17.  Growth  Rate  of  Vertically  Integrated  Global  RMS  Divergence  for  High- 
Resolution  Forecasts  Using  Different  Methods  of  Initialization,  Beginning  From 
00Z  17  January  1978  (FGGE  III-A,  Hough  Analysis) 

data  (see  Figures  18a  and  20a).  The  typical  rms  500 mb  height  changes  resulting 
from  nonlinear  normal  mode  initialization  were  12  m  for  the  Hough  analysis  and 
25m  for  the  01,  while  the  typical  rms  changes  in  the  average  1000-200mb  heights 
were  15m  and  28m  respectively.  Figure  17  shows  an  example  of  the  vertically 
integrated  global  rms  divergence  computed  from  spectral  coefficients  at  the  model 
a  -layers  as  a  function  of  forecast  time  for  the  HO-NO,  HO-DD,  and  HO-NM  fore¬ 
casts  beginning  from  the  Hough  analysis  of  00Z  17  January  1978.  TheOI  example 
is  shown  in  Figure  19  for  forecasts  beginning  from  12Z  17  February  1979. 

Returning  to  Figure  17,  it  can  be  seen  that  the  Hough  analysis  is  practically 
nondivergent.  The  forecast  model  generates  its  own  divergence  field  with  a  spin- 
up  time  of  roughly  36  hours.  The  HO-NO  and  HO-NM  runs  are  quite  close  while 
the  HO-DD  forecast  develops  a  somewhat  larger  divergence  field  after  24  hours. 
This  is  caused  by  the  difference  in  the  form  of  divergence  damping  between  the  NO, 
NM  runs  and  the  DD  run.  In  the  NO,  NM  cases,  the  divergence  is  subjected  to  a 
highly  scale-selective  V ^  diffusion  with  e-folding  times  ranging  from  261  years  for 
the  largest  scales  (n  =  l)to  41  minutes  for  the  smallest  resolved  scales  (n  =  60). 
For  DD,  the  damping  term  is  independent  of  scale,  with  an  e-folding  time  of  55 
hours  (for  forecast  times  beyond  24  hours).  The  net  result  is  that  the  larger 
scales  (for  n  <  20)  are  more  strongly  damped  by  DD,  while  the  smaller  scales 

A 

(n  >.20)  are  more  strongly  damped  by  the  V  diffusion  operator.  Thus,  the  larger 
value  of  divergence  after  day  one  in  the  HO-DD  forecast  results  from  the  increased 
short-wave  activity.  Figure  18  shows  the  vertical  profiles  of  the  global  rms 
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Figure  19.  Growth  Rate  of  Vertically  Integrated  Global  RMS  Divergence  for 
High-Resolution  Forecasts  Using  Different  Methods  of  Initialization,  Beginning 
From  12Z  17  February  1979  (FGGE  III-A,  Optimal  Interpolation) 


divergence  for  the  forecasts  (at  the  model  ff-layers)  at  (a)  the  initial  time,  (b) 

12  hours,  (c)  24  hours,  and  (d)  48  hours.  An  interesting  feature  of  the  initial 
data  is  that  the  nonlinear  normal  mode  initialization  procedure  modestly  increases 
the  divergence  in  contrast  to  the  large  decrease  found  in  the  01  data  after  the 
initialization  procedure  is  applied  (see  Figures  20a  and  21  and  discussion  below). 
Throughout  the  48 -hour  period,  the  HO -NO  and  HO-NM  profiles  are  quite  close 
to  one  another.  The  divergence  spin-up  is  most  noticeable  in  the  HO-DD  profiles, 
which  exceed  the  others  at  24  and  48  hours  for  reasons  discussed  above. 

In  view  of  the  small  differences  between  the  HO-NO  and  HO-NM  forecasts,  it 

was  decided  to  carry  out  similar  tests  with  1979  FGGE  III-A  data  which  was  pro- 

12 

duced  by  an  OI  scheme.  This  analysis  scheme  does  not  contain  any  implicit 
linear  normal  mode  initialization.  Thus,  it  was  felt  that  these  results  would  pro¬ 
vide  a  better  measure  of  the  impact  of  initialization.  Figure  19  shows  the  verti¬ 
cally  integrated  global  rms  divergence  as  a  function  of  forecast  time  for  the 
OI-NO,  OI-NM,  and  OI-DD  forecasts  starting  from  12Z  17  February  1979.  In 
contrast  to  the  Hough  case,  it  can  immediately  be  seen  that  the  OI  analyzed  diver¬ 
gence  is  much  larger  and  that  the  nonlinear  normal  mode  initialization  procedure 
reduces  the  global  value  by  more  than  40  percent.  This  is  also  quite  apparent  in 
Figure  20a.  Another  interesting  feature  here  is  that  the  OI-NM  forecast  spins 
up  slightly  from  the  initial  value  of  divergence,  while  in  the  OI-NO  forecast,  the 
model  damps  the  unusually  large  initial  value.  As  in  the  Hough  case,  the  OI-NO 
and  OI-NM  runs  are  quite  close,  while  the  OI-DD  run  produces  a  larger  diver- 
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gence  field  beyond  day  one. 

Figure  20  shows  the  vertical  profiles  of  the  global  rms  divergence  (at  the 
model  a  -layers)  for  the  OI  example  for  (a)  the  initial  time,  (b)  2  hours,  (c)  6 
hours,  (d)  12  hours,  (e)  24  hours,  (f)  48  hours,  and  (g)  72  hours.  The  most 
striking  feature  here  is  the  extremely  large  value  in  the  top  model  layer  and  the 
change  in  both  the  magnitude  and  shape  of  the  profile  after  the  nonlinear  normal 
mode  initialization  is  applied.  To  be  sure  that  this  was  not  due  to  a  data  error  in 
this  particular  data  set  and/or  a  problem  with  the  vertical  interpolation  in  the 
preprocessor,  the  global  rms  divergence  was  computed  directly  on  the  mandatory 
pressure  levels  for  12Z  18  February  1979.  This  profile  is  labelled  NO  in  Figure 
21a.  Once  again,  the  extremely  large  divergence  value  at  50  mb  can  be  seen.  The 
NM  profile  in  this  figure  is  the  result  of  the  postprocessing  (on  to  the  mandatory 
pressure  levels)  of  the  initialized  spectral  coefficients.  Figure  19b  shows  the 
vertical  profiles  of  the  model  <7-layers  and  appears  quite  similar  to  Figure  20a. 
Similar  profiles  were  also  found  in  January  1979  data. 

Returning  to  Figures  20b  -  20g  it  can  be  seen  that  the  largest  differences 
between  the  OI-NO  and  OI-NM  profiles  occur  prior  to  24  hours.  Clearly,  beyond 
24  hours,  the  combined  effects  of  the  semi-implicit  time  scheme  and  the  subgrid 
scale  diffusion  are  sufficient  to  control  the  spurious  gravity  waves  even  in  the 
no-initialization  case.  As  in  the  Hough  case,  the  OI-DD  forecast  here  also  pro¬ 
duces  values  of  divergence  that  are  larger  than  those  in  the  OI-NO  and  OI-NM 
forecasts. 

Finally.  Figure  22  shows  northern  hemisphere  plots  of  the  differences  between 
the  OI-NM  and  OI-NO  velocity  potential  at  O  =  .  5  at  (a)  6  hours,  and  (b)  48  hours 
for  the  forecasts  from  12Z  17  February  1979.  As  might  be  expected,  the  6-hour 
field  shows  quite  a  bit  of  activity,  while  the  48-hour  field  reflects  the  small 
difference  between  the  OI-NM  and  OI-NO  forecasts  at  this  time.  The  correspon¬ 
ding  500 mb  rms  height  differences  between  the  two  forecasts  are  9m  and  6m  at 
6  hours  and  48  hours  respectively.  Based  on  these  sample  forecasts,  it  is  seen 
that  the  major  impact  of  initialization  is  found  in  the  very  short  range  (less  than 
24  hours)  forecasts.  The  main  role  of  these  short  range  forecasts  (6-12  hours) 
is  to  provide  a  first-guess  field  for  the  analysis/data  assimilation  cycle.  As 
such,  the  nonlinear  normal  initialization  procedure  should  be  viewed  as  part  of 
and  necessary  for  the  analysis  scheme.  While  it  may  help  the  forecast,  it  does 
not  seem  to  play  a  crucial  role  in  shaping  the  forecasts  in  the  range  of  24-96  hours. 
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Figure  21.  Vertical  Profiles  of  RMS  Divergence  for  High-Resolution  Data  on 
12Z  18  February  1979:  (a)  Mandatory  Pressure  Levels,  and  (b)  Model  O  Layers 


9.  CONCH  SION 

The  flexible  resolution/truncation  version  of  the  AFGL  global  spectral  model 
as  adapted  to  the  CRAY-1  has  been  described.  The  numerical  formulation  and 
physical  parameterizations  are  the  same  as  in  the  baseline  model  described  by 
Brenner  et  al.  Significant  restructuring  of  the  code  was  necessary  for  efficient 
use  of  the  CRAY'S  vector-processing  capabilities.  This  resulted  in  a  substantial 
reduction  of  execution  time  for  the  forecast  model. 

To  test  the  performance  of  the  model,  a  series  of  low-resolution  (15R6)  and 
high-resolution  (30R12)  forecasts  have  been  run  and  compared.  The  12  layers  of 
the  high -re solution  model  are  the  same  as  in  the  NMC  operational  model.  The 


initial  and  verifying  data  consisted  of  FGGE  III-A  analyses  (Hough)  from  January 
and  July  1978.  In  general,  increasing  the  resolution  resulted  in  reduced  forecast 
errors  in  the  24-96-hour  range.  By  96  hours,  the  low-resolution  model  loses  its 
skill  in  forecasting  the  height  fields.  The  high-resolution  forecasts  are  still 
useful  at  96  hours  and  probably  retain  some  skill  for  forecasts  up  to  six  or  seven 
days.  The  major  weakness  is  in  the  humidity  forecasts,  which  show  very  minimal 
skill.  This  characteristic  seems  to  be  rather  insensitive  to  increases  in  resolu¬ 
tion.  This  is  caused  partly  by  the  poor  quality  of  the  humidity  fields  produced  by 
the  Hough  analysis.  Other  possible  sources  of  this  problem  are  the  moisture 
physics  schemes  and  the  vertical  advection  scheme  as  applied  to  specific  humidity. 
Both  of  these  possibilities  are  currently  being  studied. 

Two  other  topics  studied  with  the  high -re solution  model  were  the  representa¬ 
tion  of  topography  and  the  effects  of  initialization.  In  the  former,  we  found  that 
the  most  appropriate  terrain  field  for  the  forecast  model  is  the  so-called  smoothed 
topography.  It  is  produced  by  taking  the  gridded  (2.5®  x  2.5°  )  topography,  pas¬ 
sing  it  through  a  nine-point  smoother  twice,  interpolating  to  the  model’s  Gaussian 
grid,  and  spectrally  truncating  at  the  model's  resolution.  The  reasons  for  this 
choice  are  to  provide  consistency  between  the  pre-  and  postprocessing  procedures 
and  to  have  a  uniquely  defined  and  reproducible  terrain  field. 

In  assessing  the  impact  of  initialization,  it  was  found  that  beyond  24  hours, 
the  forecast  is  not  greatly  affected  by  filtering  or  initialization  of  the  initial  data. 
The  damping  characteristics  of  the  semi-implicit  time  scheme  and  subgrid  scale 
diffusion  can  control  any  undesirable  gravity  waves.  The  effects  of  the  normal 
mode  initialization  procedure  are  felt  mostly  in  the  very  short  range  (less  than 
24  hours)  forecasts,  and  it  is  therefore  needed  primarily  to  provide  smooth  first- 
guess  fields  (6-  or  12-hour  forecasts)  for  the  analysis/data  assimilation  cycle. 

Finally,  the  future  research  efforts  in  global  modeling  will  be  directed 
towards  improving  the  accuracy  and  efficiency  of  the  forecasts  with  special  emph¬ 
asis  on  improving  the  humidity  forecasts  and  implementing  cloud  forecasting 
routines.  Much  of  the  work  over  the  next  two  years  will  be  incorporating  and 
testing  the  newly  developed  physical  parameterization  packages,  which  include 
boundary  layer,  cumulus  parameterization,  and  radiation  schemes.  In  addition, 
more  efficient  transform  routines  will  be  incorporated  as  well  as  alternative 
numerical  treatment  of  the  moisture  equation.  A  more  accurate  humidity  analysis 
scheme  will  also  be  completed  within  the  next  year. 
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List  of  Symbols 


radius  of  the  earth 
gravity 

A 

indices 

vertical  unit  vector 

order  of  spherical  harmonic  (superscript) 

degree  of  spherical  harmonic  (subscript) 

pressure 

surface  pressure 

In  p* 

time 

eastward  and  northward  components  of  velocity 
nonlinear  advection  terms 
specific  heat  at  constant  pressure 
horizontal  divergence 
diabatic  heating 

meridional  derivative  of  the  associated  Legendre  function 
number  of  vertical  layers 
truncation  value  of  m 
truncation  value  of  n 

associated  Legendre  function  of  order  m  and  degree  n 

specific  humidity 

gas  constant  for  dry  air 


temperature 

basic  state  temperature 

pseudo -velocity  =  u  cos  <t>,  v  cos  <t> 

Gaussian  weight 

relative  vorticity 

absolute  vorticity 

potential  temperature 

R/Cp 

longitude 

sin 

diffusion  coefficient 
vertical  coordinate 

d  a 

vertical  velocity  = 
latitude 

velocity  potential 

streamfunction 

geopotential 

surface  geopotential 

angular  velocity  of  the  earth's  rotation 

spherical  harmonic  coefficient  of  order  m  and  degree  n 

vertical  integral 

layer  value 

level  (interface)  value 

vector 


